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Chapter 1 

Introduction 



Devices that use or produce light play an important role in modern life. Among the 
numerous daily applications of light are displays, telecommunication, data storage, 
sensors. In industry, medicine and scientific research, optical techniques are also 
absolutely indispensable. Light is used for imaging and microscopy, but also for de- 
tecting and treating diseases, analyzing chemical compounds and investigating living 
cells on a molecular level. 

In a transparent medium, like glass or air, light propagates along a straight line. 
However, as we all know from daily experience, it is impossible to see through, for 
instance, white paint or the shell of an egg. Such opaque materials have a micro- 
scopic structure that makes it impossible for light to go straight through. Figure. 1.1a 
shows schematically what happens when a beam of light impinges on a white opaque 
object: collisions with tiny particles causes the light to spread out and lose all direc- 
tionality. This process is called diffusion of light; it is comparable to the irreversible 
diffusion process that makes a drop of ink in a glass of water spread out evenly. 

Scattering and diffusion of light are huge limitations for optical imaging, but also 
severely hinder telecommunication, spectroscopy, and other optical techniques. In 
the last few decades, a tremendous effort was put in developing imaging methods 
that work in strongly scattering media.[l] That research has brought forward impor- 
tant new imaging methods like optical coherence tomography [2], diffusion tomogra- 
phy[3] and laser speckle velocimetry[4] . 

In this thesis, a new approach is taken to tackle the problem of scattering. We de- 
velop a wavefront shaping technique to steer light through opaque objects. When we 
shape the wavefront so that it exactly matches the scattering properties of the object, 
the object focuses light to a point (see Fig. 1.1b). The term 'opaque lens' was intro- 





Figure 1.1: Principle of opaque lensing. a) A plane wave impinges on an opaque scattering object. 
In the object, light performs a random walk with a typical distance given by the mean free 
path for light (£). The little light that makes it trough Is scattered in all directions, b) The 
incident wave is shaped to match the scattering In the object. The opaque object focuses 
the shaped wave to any desired point, thereby acting as an 'opaque lens'. 
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Figure 1.2: Interference in an opaque lens, a) Transmitted intensity of an unshaped incident beam. 
Scattered ligfit forms a random interference pattern known as laser speckle. Inset (phasor 
plot), at each point many waves interfere randomly, resulting in a low overall intensity, b) 
Transmitted intensity of a wavefront tfiat is optimized for focusing at a single point. The 
intensity in the focus is approximately 1000 times as high as the average intensity in a. 
Inset, in the focus all waves are in phase. 



duced in a news item[5] about our research to describe this focusing behavior. Using 
our wavefront shaping method, we steered light through opaque objects, focused it 
inside and even projected simple images through the objects. 

In this chapter, we first explain in general terms how we shape the wavefront to 
control the propagation of light. Then, we relate our work to other experimental 
fields. In Section 1.3, we give a brief introduction to the concepts and tools that are 
used for analyzing light propagation in opaque objects. We end this introductory 
chapter by giving an overview of this thesis. 



1.1 Opaque lenses 

The particles in a strongly scattering medium are smaller than the wavelength of light. 
Therefore, the wave nature of light needs to be taken into account and light propaga- 
tion cannot be described in terms of light rays. Diffusion of waves is fundamentally 
different from diffusion of particles since waves exhibit interference. We first briefly 
explain how a wave propagates through a disordered medium. Then, we show how 
interference was used in our experiments to make an 'opaque lens' focus light. 

Wave propagation in a disordered medium can be made insightful with the help 
of the Huygens principle [6]. When an incident beam of light hits a small particle 
in the object, part of the light is scattered and forms a spherical wave moving away 
from the particle. In turn, this spherical wave hits other particles, giving rise to more 
and more waves. Light propagation in a disordered scattering medium is extremely 
complex; light is typically scattered hundreds or thousands of times before it reaches 
the other side of the sample. Figure 1.2a shows the transmitted intensity for a sample 
that is illuminated with laser. The complicated random pattern is the result of the 
interference of very many different waves; this pattern is known as laser speckle. 

We now illuminate the opaque object with thousands of light beams, instead of 
one. Each of this beams forms a different random speckle pattern when it is scattered 
by the object. The total field at a given point behind the sample is the sum of the 
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speckle patterns of all incident beams. Since the object is disordered, each beam con- 
tributes to the field with a random phase. Therefore, the contributions from different 
beams mostly cancel each other (see the inset in Fig. 1.2a for a graphical representa- 
tion of this statement). 

Next, a spatial light modulator is used to delay each of the beams with respect to 
the others and thereby shape the wavefront of the incident light. The modulator is 
a liquid crystal on silicon display (LCoS). These displays have been developed in the 
last decennium for use in commercial video projectors. A computer controlled al- 
gorithm optimizes the intensity in a single target point. It does this by changing the 
phase of the beams one by one, until their speckle patterns are all in phase in the tar- 
get. In Fig. 1.2b we show the transmission after a successful optimization. All contri- 
butions interfere constructively in the target and the intensity increases dramatically 
(by a factor of thousand in this case). The opaque object now sharply focuses light 
to a single point. The shaped wavefront uniquely matches the scattering object like 
a key matches a lock. When the microscopic scatterers in the object have a different 
position or orientation, a completely different wavefront is required. Therefore, this 
method only works when the scatterers do not move during the optimization. 



1.2 Relation with earlier work 

Our work is inspired by an article of Dorokhov[8] about electron scattering in a metal 
wire. Dorokhov predicted that there exist electron wave functions that are fully trans- 
mitted through the wire, regardless of the length of the wire. Scattering of light in a 
disordered medium is, in many ways, analogous to scattering of electrons in a wire. 
The optical equivalent of this prediction is that it should be possible to construct a 
shaped wave of which all 100% of the intensity is transmitted through an otherwise 
opaque object. We wanted to construct this wavefront. The work on opaque lenses 
was initially performed as a first step towards achieving full transmission. Because 
opaque lenses turned out to be a fascinating research subject on their own, we first 
explored that field for a while. Finally, we we used our optical analogue of a disor- 
dered electronic system to confirm Dorokhov's hypothesis. 

In our opaque lens experiments, we 'borrowed' a lot of ideas from time-reversal 
experiments with ultrasound and microwaves. In pioneering work by Fink et al. (see 
Ref. (author?) 9 for a review) a short pulse impinged on a strongly scattering system. 
The amplitude of the scattered wave was recorded and played back reversed in time. 
Thanks to time reversal symmetry, the wave focused back to the original source. In 
a series of beautiful experiments, it was shown that time-reversal can be used to fo- 
cus light through a disordered medium[10], break the diffraction limit[ll, 12], and 
improve communication by using disorder[12]. 

In this thesis we show experiments that use our wavefront shaping method to ob- 
tain similar result with light. Our approach has two fundamental advantages over 
time-reversal methods. First of all, no time reversal symmetry is required. And sec- 
ondly, for time reversal one needs to have a source at the target focus. With our 
method, it is sufficient to only have a detector in the target focus. This difference 
allows us to focus waves inside a scattering medium, as is demonstrated in Chapter 5. 

An optical method that is related to time-reversal is optical phase conjugation[13]. 
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Using a non-linear crystal and high power lasers, it is possible to reverse both the 
direction and the phase of a speckle pattern to project the light back to its original 
source. Like with time -reversal, this method requires a source at the point where we 
want the waves to focus. 

Our experimental apparatus is similar to the setups used in adaptive optics [14]. 
Adaptive optics is a technique for correcting aberrations in lenses or other transpar- 
ent media, such as a turbulent atmosphere [15] or the human eye [16] . Adaptive optics 
works for situations where light propagates along rays. Diffraction and interference 
effects are not taken into account. In opaque scattering media, however, there are 
no rays of light that one can steer; light propagation is dominated by diffraction and 
interference. Therefore, our method uses interference instead of ray optics to steer 
light. 



1.3 Mathematical tools for analyzing complex system 

The propagation of light is described perfectly well by Maxwell's equations, so why 
would we need anything else? The problem is that these equations are so hard to 
solve that exact results can only be found for very simple geometries. Even for a sim- 
ple sphere geometry the result is not a closed expression, but a complicated sum of 
Bessel functions. [17] Obviously, this approach cannot be scaled to a system contain- 
ing billions of spheres, let alone irregularly shaped grains. Even worse, we do not even 
know the exact positions and orientations of the scatterers in a sample to begin with. 
Although this seems to be a hopeless situation, it is possible to capture the overall 
characteristics of the system using statistical and physical tools. In this section, we in- 
troduce the most important statistical tools that are used in this thesis: ensemble av- 
eraging, correlation functions and probability density functions. Then, we highlight 
the powerful physical concept of diffusion. These tools are applicable to all complex 
systems, whether it is light propagation in a disordered medium, or the dynamics of 
a complex biological, chemical or economical system. 



1.3.1 Ensemble averaging 

Every sample scatters light in a unique way. Even if all macroscopic properties (layer 
composition and thickness, porosity, etcetera) are the same, the microscopic struc- 
ture of two samples will be completely different. Therefore, it is impossible to predict 
the exact scattering properties of a specific sample. 

Instead, one calculates average quantities for a whole ensemble of samples. For 
example, we could calculate the optical field averaged over all conceivable samples 
that consist of a lOjjim-thick layer of zinc oxide pigment. We write this quantity as 
(E), where E is the field and (•) denotes averaging over the ensemble of all possible 
samples of a given type. 

In an experiment it is, most of the time, not needed to average over an ensemble 
of samples. Instead, one can average over the response of different portions of the 
sample. When both averaging methods are equivalent, the system is called ergodic. 
We assume ergodicity for all our samples. 
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1.3.2 Correlation functions 

The majority of the research in random scattering is concerned with calculating or 
measuring correlations of some kind. A correlation function relates the value of a 
quantity at one coordinate to the value of that quantity at a different coordinate^. For 
instance, the position correlation function of the electrical field is defined as 

C£(ri,r2)=(£*(ri)£(r2)), (1.1) 

where * denotes the complex conjugate. When two quantities are statistically inde- 
pendent, they can be averaged separately. For instance, when two points ri and r2 
are so far apart that the fields at both points are uncorrelated we can write (note the 
essential difference in the placement of the brackets) 

CB(ri,r2)=(£*(ri))(£(r2)) for ri far from rz. (1.2) 

Since E oscillates rapidly around 0, it quickly averages out. Inside a diffusive medium 
{E) ^ 0, and Eq. (1.2) vanishes. However, in general the decomposition in Eq. (1.2) 
cannot be made. Especially when ri — r2, the correlation function cannot be sepa- 
rated 

CB(r,r)= (£*(r)£(r)) / (£*(r)) (£(r)). (1.3) 

We adopt the convention that the intensity I (unit Wm^^) jg defined as / = \E\^ (see 
e.g. [19]). Using this convention, C£(r,r) equals the average intensity. Since the in- 
tensity is always positive, its average will not vanish and Csir, r) is positive. 

1.3.3 Probability density functions 

The probability density function (pdf ) gives the probability that a variable has a cer- 
tain value. An important pdf is that of the intensity / of a speckle [20] 

Eq. (1.4) tells us that the most likely intensity in a speckle pattern is zero. The chance 
of finding bright speckles decreases exponentially with intensity of that speckle. A 
typical speckle pattern with this distribution is shown in Fig. 1.2a. 

A different pdf that is used extensively in this thesis describes the field of a speckle. 
The joint probability density of the real part of the field {Ey) and the imaginary part 
of the field (£,) is given by a so-called circular Gaussian distribution 

Mi^,.E.)=— exp(-^^j. (1.5) 

The Gaussian distribution is very common and always arises when many uncorre- 
lated random variables with a finite variance are added. This important statistical 
fact is known as the Central Limit Theorem. In the case of a speckle pattern, the field 
in a single speckle is the sum of the contributions from a large number of light paths. 
When these light paths are independent, the field has a Gaussian distribution^. 



' Here 'coordinate' can denote position, angle, frequency or any other independent variable. It is also very 

common to have correlation functions involving four or more coordinates. 
^The assumption of independent paths is not completely true. Tiny deviations from Gaussian statistics 
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1.3.4 The diffusion equation 

The average propagation of intensity in a disordered medium can be described very 
well with a diffusion equation. For continuous wave illumination, the stationary dif- 
fusion equation[22, 23] applies 

-DV2/(r;f) = S(r;r), (1.6) 

where D is the diffusion constant for light, S is the source of diffuse light, and V^ is 
the Laplace operator. 

We now solve the diffusion equation to find the intensity distribution in a sample. 
AU our samples are effectively infinitely large in the x and y directions and have a 
finite thickness L in the z direction. For such a geometry, it is convenient to work in 
Fourier transformed coordinates qi = {qx,%) for the traversal coordinates x and y. 
In these coordinates, Eq. (1.6) transforms to 

^l/(q.,.) ^ = ^. (1-7) 

with (7i = |qi|. We use the Dirichlet boundary conditions 7(— z^i) — 0, and I{L+Ze2) — 
to describe the interfaces of the sample. These boundary conditions give a much 
simpler and more insightful result than the slightly more accurate mixed boundary 
conditions. [24] Here, L is the thickness of the sample and z^i and Z(,2 are the so called 
extrapolation lengths that account for reflection at the front and back surface of the 
medium, respectively. These extrapolation lengths depend on the effective refractive 
index of the sample and the refractive index of its surroundings [26, 27, 28]. When a 
slab of diffusive material is illuminated by an external source, the incident light can 
be described by a diffuse source at a depth of one mean free path £.[29] Then, the 
solution to Eq. (1.7) is (see e.g. Ref [30]) 

^^*1^^ D^isinhC^iL,) ^-^0 

/(qi,z)=-{ , (1.8) 

L^^l^^ D^isinh(^iL,) ^<^o 

with Le^ L + Zei + Ze2- A numerical Fourier transform gives I in real space coordi- 
nates. The intensity distribution in a diffusive slab is shown in Fig. 1.3. The inten- 
sity is maximal close to the source and spreads out over the medium. Far away from 
the source, the intensity decays exponentially with distance. From an expansion of 
Eq. (1.8) in small variable q±, we find that the intensity decreases exponentially with 
a decay length of Lei '/&. 



1.4 Outline of this thesis 

This thesis describes experiments and theory on controlling the propagation of light 
in opaque objects. The theoretical framework was developed in parallel with the ex- 
periments. Therefore, most of the theory is presented together with the experimental 



have been observed experimentally. [21] In Chapter 9 we describe an experiment in which we observed 
very large effects of correlations between paths. 
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Figure 1.3: Average distribution of diffuse lighit in a slab (shaded rectangle of tfiickness 30|im) of 
strongly scattering material, a) Contour plot of the energy density distribution inside the 
sample. The source is placed at a depth of z = £ = 0.72 |im and x = y = 0. The numbers 
on the equal-density curves indicate log'°(/). b) Logarithmic plot of the intensity profile at 
the bacl< of the sample (at z = 30 |im). 



results. The simple notation that was used in the first experiments did not suffice to 
describe the more advanced experiments and a more powerful notation was devel- 
oped along with the theory. Therefore there are slight differences in the notation in 
the different chapters. Care has bee taken that all chapters are self-contained and can 
be read and understood separately. 

In Chapter 2 we describe the experimental apparatus, the samples and the con- 
trol program that were used in the experiments. Special attention is given to a novel 
wavefront modulation method that we developed. 

In Chapter 3 the first experimental results of opaque lensing are presented. We also 
explain the optimization algorithm and calculate the maximally achievable intensity 
of the focus. The focusing resolution of opaque lenses is studied quantitatively in 
Chapter 4. It appears that an opaque lens focuses light as sharply as the best possible 
lens. In Chapter 5 the concept of opaque lenses is extended to focus light inside an 
opaque object. In Chapter 6 we demonstrate that disorder can improve optical com- 
munication. In Chapter 7 dynamic algorithms for use with non-stationary samples 
are investigated. 

In Chapter 8 a transport theory for optimized wavefronts is developed. We show 
that wavefront shaping significantly increases the total transmission through opaque 
objects. The experimental observation of this effect is presented in Chapter 9. 

Finally, in Chapter 10, we summarize our findings and give an outlook of the many 
possible applications of our work. 
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Chapter 2 

Experimental apparatus 



In this chapter we discuss the experimental apparatus that we built to control the 
propagation of light in disordered media. We discuss the considerations that played 
a role in designing the experiment and give special attention to a new wavefront shap- 
ing method that we developed for our experiments. The setup was modified for each 
of the individual experiments that are described in this thesis. Here we discuss the 
common elements of the design and explain what specific modifications were made. 

A general diagram of the experiments is shown in Fig. 2.1. A strongly scattering 
sample is illuminated by a wavefront synthesizer that is able to construct a spatially 
modulated beam. The light that is scattered by the sample is collected by a detec- 
tion system. This detection system provides feedback to a computer algorithm that 
controls the wavefront synthesizer. 

The most complex part of the apparatus is the wavefront synthesizer. It contains a 
liquid crystal display (LCD) that spatially modulates incident light. In Section 2.1 we 
analyze the characteristics of the LCD and explain how it was used as a phase modu- 
lator. We also developed a new method for using the LCD to spatially modulate both 
amplitude and phase of the light. This method was published in Ref. (author?) 1. In 
Section 2.2 the detection system is described. Here we also discuss how the detection 
is synchronized with the wavefront synthesizer. We comment on the overall stability 
issues of the setup in Section 2.3. Then, in Section 2.4, we describe which types of 
samples were used and introduce a new fabrication method that we developed for 
making strongly scattering samples. The structure of the control program that coor- 
dinates the experiment is discussed in Section 2.5. Finally, in Section 2.6 we give an 
outlook of what can be expected in future experiments. 



wavefront 
synthesizer 



control 
sample detection program 




Figure 2.1 : Block scheme of the experiment. A wavefront synthesizer generates a shaped monochro- 
matic wavefront. The light is scattered by a strongly scattering sample. A detector defines 
a target position for the optimization procedure and provides a feedback signal. A computer 
analyzes the signal and reprograms the phase modulator until the light optimally focuses 
on the target. 
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Figure 2.2: Operating principle of a transmissive TN LCD. a) No voltage is applied and the rod-like 
liquid crystal molecules are ordered in a helix. The polarization of incident light follows the 
twist, b) A voltage is applied between the two electrodes. The rods orient in the direction 
of the field and the twist disappears. The polarization of the incident light is not rotated. 

2.1 Wavefront synthesizer 

Computer controlled wavefront shaping is a very versatile technique that is used in 
many fields of optics. In adaptive optics [2], for example, deformable mirrors or other 
spatial light modulators are used to correct aberrations in a variety of optical sys- 
tems. Another area that relies on computer controlled light modulators is that of dig- 
ital holography. Topics in digital holography include holographic data storage [3], 3D 
display technology [4], and holographic image processing[5]. 

Liquid crystal displays are among the most popular types of light modulators be- 
cause of their high optical efficiency, the high number of degrees of freedom and their 
wide availability. In our experiment, light is modulated with a twisted nematic (TN) 
liquid crystal display^. This LCD can modulate light at a refresh rate of 60 or 72 Hz at 
a resolution of 1024 x 768 pixels. In Section 2. 1 . 1 we discuss the operating principle of 
a TN LCD. LCDs are designed to modulate light intensity. To use the LCD as a phase 
modulator, a thorough characterization of the display is required. This characteriza- 
tion procedure is described in Section 2.1.2. In Section 2.1.3, we describe a common 
method to use TN LCDs as phase modulators. This method was used to perform the 
experiments that are described in Chapters 3 and 7. For the other experiments in this 
thesis, we developed a new modulation method that has several advantages over ex- 
isting methods. This novel modulation technique is introduced in Section 2.1.4 and 
experimentally demonstrated in Section 2.1.5. The switching behavior of the display 
is described in Section 2.1.6. 



2.1 .1 Principle of a twisted nematic liquid crystal display 

The operating principle of LCDs is based on the birefringence of rod-like organic 
molecules. In the twisted nematic phase, these rods are ordered in a helbc as is shown 
in Fig 2.2a. A transmissive LCD is typically sandwiched between two crossed polar- 



'We used two Holoeye LC-R 2500 systems. The LCDs have a size of 19.6min x 14.6mm. One system was 
customized by the manufacturer to control two LCDs with one control box. That system was used for 
the experiments in Chapter 9, where we needed control over both polarizations. 
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izers. The 90° helical structure of the liquid crystal rotates the angle of polarization 
of the incident light so that the light passes the second polarizer. When a voltage is 
applied over the liquid crystal cell, the molecules align with the electrical field, as is 
shown in Fig 2.2b. Now the optical axis of the molecules is parallel to the direction of 
light propagation and the polarization of the incident light is not rotated. As a result, 
the light is blocked by the second polarizer and the pbcel is dark. 

The LCDs in our experiments are TN liquid crystal on silicon (LCoS) displays. These 
reflective displays are designed to be used with a polarizing beam splitter or with 
crossed polarizers in an off axis configuration. The thickness of the liquid crystal 
layer and the twist angle are carefully engineered to optimize brightness, contrast 
and response time for projecting video. 

The operation principle of a reflective LCD is more complicated than that of a 
transmissive LCD (see e.g. Ref. (author?) 6). When no voltage is applied, the polariza- 
tion of the linearly polarized incident light follows the helix. At the back surface of the 
LCD the light is reflected and, on its way back it is rotated back to its original polar- 
ization. When a voltage is applied, the helix is distorted (see Fig. 2.2b). Now, linearly 
polarized light cannot completely follow the twist anymore and becomes elliptically 
polarized. In the 'on' state, the light is exactly circularly polarized at the back surface 
of the LCD. Reflecting circularly polarized light changes its handedness. Therefore, 
on its way back, the polarization is not rotated back to its original polarization, but 
to the orthogonal polarization state. In the crossed polarizer configuration, the pixel 
now is reflecting. 

2.1.2 Liquid crystal display characterization 

The optical characteristics of a single pbcel of the LCD can be described by its Jones 
matrbc /. The Jones matrix relates the horizontal and vertical polarization compo- 
nents of the incident field (denoted as, respectively EJJ and £'") to the components 
of the outgoing field 
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Elliptic and circular polarization are described by complex values of Eh and Ey. In 
general, the elements of the Jones matrix are also complex numbers. A pixel of the 
LCD is fully characterized by measuring J[V) for all voltage settings V. 

The setup in Fig. 2.3a was used to measure the Jones matrix of the LCD. A laser 
illuminates a disk with a diameter of approximately 3 mm in the center of the LCD. 
With two waveplates, any polarization of incident light can be generated. A third 
waveplate and a polarizer are used to analyze the modulated light in any desired po- 
larization basis. The modulated light is focused on a detector. To obtain, for instance, 
the /12 component, the polarization optics were rotated so that the incident light was 
vertically polarized and the reflected horizontally polarized light was analyzed. We 
later used the simpler setup that is shown in Fig. 2.3b. That setup has no polarization 
optics that need to be rotated and, therefore, is less sensitive to alignment inaccu- 
racies. With the simplified setup only the Ju component can be measured. For the 
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Figure 2.3: Setups used for characterizing the LCD. a) Setup for full characterization. A/2, half wave- 
plate; A/4, quarter waveplate; BS, non-polarizing 50% beam splitter; P, polarizer; L, 60cm 
lens; D, detector. The light source was a 632.8 nm HeNe laser, b) Setup for partial charac- 
terization. PBS, polarizing beam splitter. Partial characterization was performed at wave- 
lengths of 632.8nm (HeNe laser) and 532nm (Nb:YAG laser). 



modulation scheme that is discussed in Section 2.1.4, such a partial characterization 
of the LCD is sufficient. 

In both setups, /12 was measured using a diffractive technique that is comparable 
to the method described in Refs. (author?) 7 and (author?) 8 but only required detec- 
tion of the 0* order diffraction mode. We assumed that all pixels of the LCD have the 
same Jones matrbc. First the absolute value of 712(1^) was obtained by varying the volt- 
age over all pixels of the LCD simultaneously. Then, we programmed the LCD with a 
binary grating with a duty cycle of 50% (a so called Ronchi grating) . The voltage over 
the notches of the grating was kept constant at \o while the voltage over the rules was 
varied. The intensity in the 0* diffraction order responds as 






= l/i2(t^o)l' + 1/12(^^)1' +2|/i2(Vo)ll/i2(V)| COS [arg/i2(V)-arg/i2(Vo)], (2.2) 



which gives us, in principle, enough information to obtain the phase of JiziV) up to 
an overall phase offset. In practice, however, the inversion of Eq. (2.2) is very sensi- 
tive to noise for some values of V. We solved this problem by measuring I{V, %) for 
different values of Vq (0%, 25%, 50%, 75% and 100% of the maximum voltage). Then 
we combined the data, using only measurements for which the inversion of Eq. (2.2) 
was stable. The result of these measurements is shown in Fig. 2.4. We find that the 
phase of the reflected light changes from tt to with increasing voltage. The ampli- 
tude increases from to a maximum at the 'on' state (around a phase of 35°) and then 
decreases slightly. 

To obtain all four elements of the Jones matrix, the measurement was repeated for 
each element of /. Moreover, the measurements were also performed in a rotated 
basis with left and right hand circularly polarized light. These extra measurements 
resolved the relative phase of the different components of the / matrix. 

All in all, measuring the full Jones matrix of an LCD is cumbersome. A complicating 
factor is the fact that all polarization states (except for horizontal or vertical polariza- 
tion) are changed by reflecting off a coated mirror or passing through a beam splitter. 
Theoretical analysis is complicated by the fact that in thin LCDs surface effects start 
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Figure 2.4: Polar plot of the modulator output amplitude witin vertical polarization in, horizontal polariza- 
tion out (/12), at a wavelength of 532 nm. The modulation voltage increases in the direction 
of the gray arrows. 



to play a role. Also, the total amount of light that was reflected by the LCD was found 
to depend on the voltage, which means that the Jones matrix is not unitary. There- 
fore, many of the theoretical models that are used to describe LCDs[9, 10, 1 1] can only 
give an approximate result. 



2.1.3 Liquid crystal phase-mostly modulation 

TN LCoS displays for intensity modulation typically achieve a maximum phase re- 
tardation of around n for horizontally or vertically polarized light. Nevertheless, by 
choosing a suitable combination of incident polarization and analyzer orientation, it 
is often possible to find an operating mode where the phase retardation is 2n while 
the amplitude modulation is relatively low. This mode of operation is called 'phase- 
mostly' modulation and has been the subject of extensive experimental and theoret- 
ical research[7, 8, 10, 11, 12, 13, 14, 15, 16, 17, 18]. 

The optimal configuration of the polarizers is different for each specific series of 
LCDs and also depends on the wavelength of the light. In general, the optimal com- 
bination of polarization states is elliptical. [17] We found the optimal polarizations for 
a wavelength of 632.8 nm with a brute force optimization method. This method used 
the measured Jones matrix of the LCD to numerically try all configurations. In the 
optimal configuration the transmission is low and the system is very sensitive to vari- 
ations in the angles of the waveplates. The measured modulation curve (see Fig. 2.5) 
has an intensity variation of21%ofthe mean value. This configuration was used suc- 
cessfully in the wavefront shaping experiments that are described in Chapters 3 and 7. 
In these experiments, a reference detector compensated for the intensity variations. 
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Figure 2.5: Polar plot of the modulator output amplitude. Dashed circle, perfect phase only modulation; 
solid curve, measured response in the configuration for optimal phase-mostly modulation 
at a wavelength of 632.8 nm. 



2.1.4 Decoupled amplitude and phase modulation^ 

In this section, we describe a novel modulation scheme that uses a TN LCD in combi- 
nation with a spatial filter to achieve individual control over the phase and the ampli- 
tude of the light. The method has four major advantages over the scheme described 
in the previous section. Firstly, it allows separate control over amplitude and phase 
of the wavefront. Secondly, the residual amplitude-phase cross-modulation that oc- 
curred with the phase-mostly modulation scheme is almost completely eliminated. 
Thirdly, the LCD operates in a simple horizonal-in vertical-out configuration, this 
saves components and characterization time, and makes the setup less sensitive to 
align. And finally, the experimental setup can be easily extended to control both po- 
larizations. 

Since the introduction of liquid crystals, many techniques have been developed to 
achieve combined amplitude and phase modulation. Examples are setups where two 
LCDs are used to compensate amplitude-phase cross modulation[19, 20], and double 
pixel setups where two pixels are combined to a macropixel[21, 22, 23]. Each of this 
techniques has its own limitations. The use of two LCDs introduces alignment and 
synchronization issues. Dual pixel schemes put specific demands on the modulation 
capacities of the LCD, such as requiring phase-only modulation[23], amplitude-only 
modulation[21] or a phase modulation range of 2?! [22]. 

We developed a method that combines four pixels into a macropixel. With this 
method, full spatial amplitude and phase modulation can be achieved with any LCD. 
The setup used for this modulation scheme is shown in Fig. 2.6. A monochromatic 
beam of light is incident normal to the SLM surface. The modulated light is reflected 
from the SLM. We choose an observation plane at which the contribution of each 



^This section and tlie following section are based on the article Spatial amplitude and phase modulation 
using commercial twisted nematic LCDs, E. G. van Putten, I. M. Vellekoop, and A. E Mosk, accepted for 
publication in Applied Optics (2008). 
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Figure 2.6: Experimental setup to decouple phase and amplitude modulation. Four neighboring pixels 
form one macropixel which can modulate amplitude and phase. In the plane of observation 
neighboring pixels are n/l out of phase (inset). A spatial filter combines the pixels into one 
macropixel. Li, L2, lenses with a focal distance of 250mm and 200mm, respectively. 



neighboring pixel is n/l out of phase, as is seen in the inset. A spatial filter removes 
all higher harmonics from the generated field, so that four neighboring pixels are 
merged into one macropixel. 

By choosing the correct combination of pixel values, any complex value of the total 
field can be synthesized. The electric field in a macropixel, Egp, can be written as the 
sum of the fields, £1, £2, E^: and £4, coming from the four different pixels. Behind the 
spatial filter, Esp is given by 



-sp ■ 



: El exp 



-in 
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+ E2 exp(/7r) + £3 exp | -in\+ £4, 



: [Eir - Ezr) + i{E3r " ^Ir) + (^i; " £3;) + iiE^ - Ezi), 



(2.3) 

(2.4) 



where the indices r and i refer to the real and the imaginary part of the field. The 
voltages on pixels 2 and 4 are chosen such that £4, — Ezt — 0, and in the same way 
the voltages on pixels 1 and 3 are chosen such that En — Est — 0. Equation 1 is now 
reduced to 

£sp = {Eir - E2r) + i[E3r " Ei,-] (2.5) 



The separate pixels are programmed so that the fields are given by 

£1 =—Asm{i^]+ iAi, 
E2 = —Acos{i^)+iA2, 
£3 —Asm{i^)+ iAi, 
£4 = Acos(j0)+ /A2, 



(2.6) 
(2.7) 
(2.8) 
(2.9) 



with A and (p, respectively, the desired amplitude and phase. The imaginary parts A 
cancel. The desired complex value is now synthesized by the real values of the field at 
the four different pixels. From the geometrical construction shown in Fig. 2.7 it can 
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Figure 2.7: Modulation amplitude response of four pixels that form a macropixel. Pixels have a njl 
phase shift with respect to each other. The four pixels synthesize any complex value: £2 
and En generate the real part of the field; lm(iJ4-£2) = 0. Ex and E-i form the imaginary 
part (not shown). 



be seen how we modulate a value on the real axis, E,^y — Ezr, by choosing £4,- — Ezt. 
It is always possible to find pixels values with exactly opposite imaginary parts and 
different real parts. The only requirement posed on the SLM is that at least one of the 
field components has to vary when the pixel voltages are changed. 

The decoupled amplitude and phase modulation was used with a A = 532 nm diode 
pumped solid state laser^ for the experiments in Chapter 5. For the experiments in 
Chapters 4, 6 and 9, a A = 632. 8nm helium neon laser was used"*. In Chapter 9 the 
wavefront synthesizer was extended to modulate two polarizations by simply replac- 
ing the beam dump (see Fig. 2.6) by a second LCD. 



2.1.5 Demonstration of amplitude and phase modulation 

We tested the new modulation technique with the same characterization method as 
was used to measure the modulation curve of the LCD (see Section 2.1.2), with the 
difference that we now use macroprxels instead of actual pixels. We programmed the 
modulator so that the macroprxels formed a Ronchi grating. The notches were set at 
half of the maximum amplitude with a phase offset of zero. The phase and amplitude 
of the rules was varied. A detector recorded the light intensity in the 0* diffraction 
order^ . These experiments were performed at a wavelength of 532 nm. 

We observed that the intensity in the O**^ diffraction order varied as the cosine of the 
phase difference between the notches and the rules of the grating, just as is expected 
from Eq. (2.2). We repeated the experiment with different amplitudes in the rules 
of the grating. All experimental result overlaps almost perfectly with the theoretical 
curves (see Fig. 2.8). From this agreement, we conclude that a full 2n phase shift is 



^Coherent Compass M315-100, max. 100 mW, intra cavity doubled, diode pumped Nb:YAG laser. 

"IDS Uniphase 1125/1? 5mW polarized HeNe laser. 

^The 0* diffraction order of the macropixel grating is the direction of the modulated light when all 
macropixels are set to the same amplitude and phase. This is not the same as the 0* diffraction order 
of the modulator, which is the direction of the light when all actual pixels are set to the same voltage. 
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Figure 2.8: Intensity in 0'*^ order grating mode as a function of tlie set phase difference A9 = dB- 6a- 
Notches of the grating are kept at amplitude A — 0.5, phase Oa — 0. Amplitude B and 
phase 0B of the rulers is varied. Soiid curves, expected intensity for perfect modulation. 
Intensities are referenced to /o = 4.56- 10^ counts/second. 




Figure 2.9: Independent phase and amplitude modulation. Curves show the measured relative ampli- 
tude A/y% as a function of the programmed phase. Relative amplitudes are set to 0.25, 
0.5, and 0.75. /o = 19.7-10=^ counts/second. 



achieved. Moreover, at n phase shift the intensity in the 0* order vanishes, which 
indicates that the field of the notches has the same magnitude and opposite sign as 
the field of the rules of the grating. 

To determine the amount of amplitude-phase cross-modulation quantitatively, we 
measured the intensity in the 0'*^ diffraction order when all macropixels are set to the 
same field and amplitude. When the phase is cycled from to 2?!, the observed inten- 
sity remains virtually equal. From the results shown in Fig. 2.9, we find that the am- 
plitude is constant to within 2.5%, which is a significant improvement over the 21% 
cross-modulation observed with the phase-mostly modulation scheme. The small 
amplitude variations are periodic with a n/2 period. This periodicity is understand- 
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able, since after 7r/2 rotation of the phase, the same pattern is programmed on the 
LCD shifted by one pixel (see Eqs.(2.6)- (2.9)). Unlike the phase-mostly modulation 
method, and 2n are equivalent, which means that there are no discontinuities at 
phase wraps. The amplitude variations are probably the result of dynamic cross-talk 
between neighboring pbcels. This effect is discussed in the following section. 

Although the amplitude of a macropbcel is not affected by its phase, neighboring 
macropbcels do effect each other. A detailed analysis of this effect is found in Ref. (au- 
thor?) 24. The effect is similar to ordinary diffraction: the beams originating from 
the two macropbcels expand through diffraction and interfere with each other. It was 
found that, typically, a completely random wavefront carries approximately 85% of 
the intensity of a plane wave. The rest of the intensity is clipped by the iris diaphragm. 
In Chapter 5 this effect was compensated for by measuring a reference intensity with 
a random wavefront. In Chapter 9 we used a more accurate method where the sample 
is translated. 

2.1.6 Transient behavior 

We now investigate how fast a pixel of the LCD can be switched. The image on the 
modulator is updated with a refresh rate of 60 frames per second. To allow for real 
time operation, we reconfigured the gamma lookup table of the LCD driver electron- 
ics. The table was configured so that a pixel value linearly corresponds to an ampli- 
tude modulation on the real or imaginary axis. Furthermore, the conversion from 
a polar representation (amplitude and phase) to the real and imaginary parts (see 
Eqs. (2.6) -(2.9)) is performed in real time by the video acceleration hardware. In the 
process of displaying a new image the following sequence of events can be identified: 

1. The control program loads new matrices for the amplitude and phase to the 
video hardware 

2. The video hardware scales and translates the matrices to screen coordinates. 
Then it performs the necessary calculations to convert amplitude and phase to 
pixel values and prepares the new image in a background buffer. 

3. The control program waits to just before the start of a new frame (the so called 
vertical retrace period) to swap the background buffer with the foreground buff- 
er. 

4. After every vertical retrace the video hardware sends the image to the light 
modulator over a digital visual interface (DVI) link. 

5. The light modulator hardware receives the image and converts pixel values to 
voltages with the use of an internal gamma lookup table 

6. The light modulator hardware drives a matrix of transistors on the back of the 
LCoS display according to a pulse width modulation (PWM) scheme. 

To measure the transient behavior of the LCD, we repeatedly switched the whole dis- 
play from the minimum to the maximum voltage and back, each time waiting 500 ms 
between switches. A camera was configured to measure the intensity in the center 
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Figure 2.10: Transient intensity response of the LCD for switchiing the amplitude from to maximum 
(solid curve) and for switching it back to (dashed curve). During the first 33 ms (2 frames) 
the display does not respond at all. The intensity is normalized to the average intensity in 
the 'on' state. 



of the modulated beam with a shutter time of 1 ms. The delay between the vertical 
retrace and the camera trigger was varied from ms to 120 ms. 

Figure 2.10 shows the measured switching response of the LCD. We recognize three 
distinct periods. During the first two frames (from to 33 ms) the image on the LCD 
does not change. We call this period the idle time Tidie. During the idle time the 
image is transferred from the computer to the modulator. In the second period the 
voltage over the pixels changes and the liquid crystal molecules reorient, which takes 
a certain time Tsettie '^ 50 ms. In the third period, the liquid crystal molecules have 
reoriented completely. However, the signal still oscillates with a period of exactly half 
a frame (8.3 ms). These oscillations are the result of how the LCD hardware drives 
the pbcels (see e.g. [25]). Each pbcel is switched on and off according to a PWM code. 
A storage capacitor at each pixel integrates the total current to achieve the desired 
average voltage over the liquid crystal. After half a frame, the controller reverses the 
voltage to avoid a DC current that would damage the liquid crystal. This switching 
scheme results in rapid oscillations in the reflected light. 

When two neighboring pixels have a different voltage, there is a field gradient at 
their border. Due to the pulse modulation scheme, the gradient will oscillate. The 
magnitude of this undesired effect depends on the PWM code of each of the pixels. 
For example, consider the field response for a varying phase and a constant ampli- 
tude of 0.25 (the smallest circle in Fig. 2.9). The response curve shows small jumps in 
the amplitude at phase values of — 15° and +15°, these jumps are repeated every 90°. 

To understand these jumps, we examine the bit code of adjacent pixels in a macro- 
pbcel. For a phase of 14°, pixel 1 and 2 have a value of 255 and 230. At a phase of 15°, 
the first pixel values has changed to 256 and the second pixel is still at 230. Although 
the change from 255 to 256 appears to be small, these numbers have a completely dif- 
ferent bit-pattern (01 II 111 11 and 100000000 respectively). Therefore, the field gra- 
dient between pixel I and pixel 2 wUl differ significantly between the two situations, 
resulting in a jump in the field response. 

The transient switching characteristics of the LCD and the temporal oscillations 
put special demands on the timing of the detection. This issue is addressed in Sec- 
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Figure 2.11 : An imaging telescope is used to project the surface of the LCD to the entrance pupil of 
the microscope objective. Crosses indicate focal points of the lenses. 



tion 2.2. 

2.1.7 Projecting the wavefront 

In most of the experiments in this thesis, the shaped wavefront is focused on the 
sample by means of a microscope objective. The surface of the LCD is imaged onto 
the entrance aperture of the microscope objective with an imaging telescope. When 
combined phase and amplitude modulation is used {see Section 2.1.4) there is a pin- 
hole in the focus of the telescope to spatially filter the generated wavefront. The tele- 
scope also serves to demagnify the beam coming from the LCD to the size of the mi- 
croscope objective's aperture. It is essential that a telescope with two positive lenses is 
used so that beams that leave the LCD at an angle are also imaged onto the aperture 
of the objective (see Fig. 2.11). When the lenses are aligned properly, the entrance 
pupil of the objective corresponds to a circular area on the LCD. Pixels outside this 
area do not contribute to the field on the sample and can, therefore, be skipped in the 
optimization procedure. 

In this configuration, a pixel on the LCD corresponds to an angle in the focal plane 
of the objective. When we group pixels on the LCD together in blocks, we reduce the 
angular resolution at the sample surface and, thereby, reduce the size of the projected 
spot. The aperture of the objective is always filled completely. Therefore, when a high 
NA objective is used, the number of segments in the incident wavefront is approx- 
imately equal to the number of mesoscopic channels on the sample surface. If, for 
example, all pixels are grouped into a single segment, light is focused to a diffraction 
limited spot encompassing exactly one mesoscopic channel. 

If the sample is not positioned in the focal plane of the objective, or when a low 
NA objective is used, the wavefront synthesizer illuminates a spot that supports more 
scattering channels than there are control segments in the incident wavefront. Thus 
the incident field cannot completely be defined by the wavefront synthesizer. It turns 
out that, in general, this limitation has little effect on how well the propagation of 
light is controlled (see for instance Fig. 3.4). 



2.2 Detection 

During the wavefront optimization procedure, a detector monitors the intensity in 
the target. Different detectors were used for this purpose. In this section we discuss 
the detectors that were used, as well as their relevant properties. 
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The first experiments were performed with a photodiode. A photodiode has an 
excellent dynamic range and a fast response. The major drawback of using a photo- 
diode is that it is very hard to select exactly a single speckle. For this reason we started 
using cameras. During optimization the camera image is integrated over a disk with a 
software defined radius and position. The diameter of the disk is chosen to be slightly 
smaller than the speckles that are visible on the initial camera image. Using a camera 
also allowed us to easily define multiple targets and to monitor the intensity in the 
background around the optimized speckle. 

We have used three different models of cameras. The first model is the Allied Vision 
Technologies Dolphin F-145B. This camera is an all purpose charge coupled device 
(CCD) camera that is connected to the computer with a IEEE 1394 link (firewire). To 
increase the dynamic range of the camera, the shutter time was varied by the con- 
trol program. In Chapter 9 the required dynamic range was so high that changing the 
shutter time was not sufficient. Instead, the computer controlled a motorized trans- 
lation stage that automatically placed a neutral density filter in front of the camera. 

We used the Dolphin camera for all experiments, except for the fluorescence mea- 
surements described in Chapter 5. For that experiment, we started with a Hama- 
matsu ORCA electron multiplying CCD (EMCCD) camera. This camera is cooled with 
a Peltier element to reduce the dark current. Also, the signal is amplified on the CCD 
chip to overcome readout noise when the signal is very weak. The Hamamatsu cam- 
era was connected to a dedicated computer with a CamLink interface. However, al- 
though the Hamamatsu camera has a high sensitivity, it was not possible to measure 
small variations in the signal intensity. It turned out that the overall background in- 
tensity of the camera image varied on a frame by frame basis. This so called base- 
line drift problem was solved by using a different EMCCD camera. The Andor Luca 
DL658M is a cooled EMCCD camera that connects to the USB 2.0 bus. It has a base- 
line clamping feature that eliminates almost all baseline drift. 

The linearity of all cameras was confirmed experimentally. We also recorded a 
background image for every experiment. For most experiments it was sufficient to 
subtract the average value of the background image from the measured signal. How- 
ever, the experiments in Chapter 9 required a measurement of the exact intensity 
distribution over the whole camera. Therefore, the full background image was sub- 
tracted pixel by pixel. Moreover, in that experiment we corrected for the approxi- 
mately 30% lower sensitivity of the camera close to the corners of the CCD chip. The 
lower sensitivity is probably the result of a minute misalignment of the microlenses 
on this chip. 



2.2.1 Timing 

In all our experiments it turned out that the speed at which the optimal wavefront 
can be constructed is the limiting factor. Therefore, we want to measure as quickly as 
possible. In Section 2.1.6 we saw that the wavefront oscillates with a period of 8.3 ms. 
To avoid noise due to aliasing, the shutter time of the cameras (TuieasJwas always set 
to a multiple of this value. 

In Section 2.1.6 we observed that it takes some time for the image on the LCD 
to start changing (Tidie) and then it takes some more time for the image to stabilize 
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Figure 2.12: Timing diagram of a series of measurements. Time (in frames of 16.7ms) is indicated 
on tine top axis. The curve in the topmost box symbolizes the switching behavior of the 
modulator. The numbered rectangles in the box below stand for the time that the camera 
is recording an image. In the lower box the relevant time intervals for synchronization are 
drawn. 



(^settle)- To maximize the number of measurements per second and to further reduce 
aliasing effects, we synchronized the camera with the wavefront synthesizer. 

A timing diagram of a series of measurements is given in Fig. 2.12. In this diagram, 
Tidie = 33 ms — 2 frames and Tsettie = 57ms — 3.4 frames. The first image is presented 
during the vertical retrace period of the modulator (frame 0). At frame 2, the image 
starts to change. Then, at frame 5, the second image is presented, although the mea- 
surement for the first image has not even started. Since it takes two more frames for 
the LCD to react, the image on the LCD is stable from frame 5 to frame 7. We trigger 
the camera after 7idie+ Tsettie — 5.4 frames. Even if there is some jitter in the timing, the 
camera has finished measuring before the image on the LCD starts changing (frame 
7). 

With this tight timing scheme, we can perform a measurement each 5 frames (in- 
stead of each 7 frames) . The idle time Tidie only affects the first measurement of a 
sequence. In the experiments we perform a sequence of 5 to 10 synchronized mea- 
surements for each segment. 



2.3 Stability 

The optimized wavefront is unique to a sample. When the sample is moved, the opti- 
mized wavefront does not fit the sample anymore. Therefore, during the course of an 
optimization, the sample has to be stable with sub -wavelength accuracy. The stabil- 
ity demands for the rest of the apparatus are high as well. In this section we discuss 
the most important stability considerations for our experiments. 

The whole apparatus was built with the requirements of interferometric stability 
in mind. We built it on an actively levelled, damped optical table^. Moreover, only 



^TMC CleanTop II 780 series, 783-655-12R with 14-426-35 support. 
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high quality opto-mechanical components were used^. The sample was mounted on 
a flexure stage to avoid the creep problems that are intrinsic to stages with ball bear- 
ings. The mechanical stability of the system was tested by tapping the optical compo- 
nents while monitoring the position of the focused light with sub-micron accuracy. If 
the focus did not return to exactly the same position after tapping a component, the 
mount of the component was replaced by a more stable one. 

The temperature and humidity of the setup were not controlled. Therefore, it is 
likely that thermal drift and hygroscopic expansion negatively affect the stability. Af- 
ter switching on the system, it takes a few hours to stabilize completely. During this 
time the laser, the cameras, the LCD, and the computer controlled translation stages 
warm up. After this warmup time the total system is, typically, stable enough to per- 
form an optimization that takes an hour. 

We found that air turbulence causes fluctuations in the feedback signal. These 
fluctuations are in the order of a few percent in the interesting frequency range of 
about 1 Hz. To reduce air turbulence, a box was built around the experiment. This box 
reduces stray light as well. As the laser is the most important source of hot, turbulent 
air, it needed to be placed outside the box. Of course, all components with fans (the 
cooled camera and the LCD control electronics) were also placed outside the box. 

Finally, there is a constraint of the stability of the wavelength of the laser. To calcu- 
late the required stability, we estimate the average path length through a sample. For 
example: a typical sample has a mean free path of the path ofl — 0.7 ^im, and a thick- 
ness of L= I5i. Then, the diffuse path length 5 is in the order of 5 — 15^i — 160 ^xm. At 
a wavelength of 532 nm, a typical path is 300 wavelengths long. When a wavelength 
change results in a ti phase shift over this path length, the wavefront optimization 
fails. In this example, the wavelength needs to be stable to within a nanometer (or ex- 
pressed in inverse centimeters, better than 31 cm i), which is absolutely no problem 
for a temperature stabilized Nb:YAG laser or a HeNe laser. 



2.4 Samples 

For our experiments we used a large diversity of strongly scattering objects. In many 
ways, our method for controlling the propagation of light in such objects does not 
depend on the optical parameters of the sample. There are, however, some practical 
limitations that restrict what samples we can use: 

Stability The optimal wavefront for focusing light through a strongly scattering ob- 
ject uniquely depends on the exact configuration of the scatterers. Therefore, 
we are limited to solid samples. We observed that after optimization the signal 
decreases with a typical timescale of approximately one hour, a value that is 
probably limited by thermal drift in the setup rather than by the samples. 

Absorption For our method to work, there has to be at least some light on the de- 
tector before the optimization procedure is started. If the absorption in the 



^Most opto-mechanical components were manufactured by Siskiyou and Thorlabs. For tlie beam ex- 
pander, the laser, the LCDs, and the sample, special stable mounts were designed. 
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sample is too strong, there will be no feedback signal to optimize. For most ex- 
periments it should be possible to use weakly absorbing materials as long as the 
initial transmission is detectable. However, for the experiment in Chapter 8, we 
expect that absorption decreases the effect that we are interested in. Therefore, 
we only used non-absorbing materials in all our experiments. 

Thickness The optimization method detects the intensity in a single scattering chan- 
nel. The number of independent channels is roughly equal to (21)^/(^/2)^, 
where L is the thickness of the sample and A is the wavelength. Moreover, the 
total transmission scales as ln/L, with £tr the transport mean free path for light 
in the medium. Therefore, the total optical power in a single speckle scales with 
L^^. For this reason, most of our samples were made relatively thin (~ 10 ^xm). 
However, successful optimizations have been performed on samples that were 
up to 1.5 mm thick (a baby tooth, see Table 3.1). 

Scattering strengtii We want to study light propagation in non-absorbing opaque 
scattering media where all transmitted light is diffuse. A scattering object is 
opaque when it is thicker than a few transport mean free paths. We require that 
L > 4£tr- In this regime, the fraction of ballistic (non-diffuse) transmission is 
less than exp(— 4) ^ 2%. To keep the number of channels as low as possible, 
it is advantageous to use thin samples with a low £tr- There are two ways to 
minimize £tr- First of all, the scale of the disorder should be comparable to the 
wavelength of light in the medium. Secondly, the index contrast should be as 
high as possible. Good candidates for making strongly scattering samples are 
pigment particles made of a high index material (Ti02, ZnO). These particles 
have a typical size of ~ 200 nm. 

Flatness and homogeneity Our method for controlling propagation of light does not 
depend on the flatness or the homogeneity of the samples. However, for sys- 
tematically analyzing the experimental results it is highly desirable to have sam- 
ples that have the same thickness over the whole sample area (< 20% varia- 
tions). Furthermore, the samples should be homogeneous in composition and 
certainly not have any holes. 

Special requirement: doping For the experiment that is described in Chapter 5 it was 
required to place fluorescent nanospheres in the scattering medium. 

Special requirement: substrate The experiment in Chapter 9 required that the sam- 
ples were on a thin glass substrate to allow two high NA microscope objectives 
to focus both on the front and the back of the sample. Working without sub- 
strate all together was not possible because the substrate provides structural 
stability to the sample. 

Special requirement: effective refractive index The experiment in Chapter 9 also re- 
quired the number of channels to be as low as possible. Therefore, we used 
thin, strongly scattering samples. Moreover, we chose to use ZnO pigment in 
an air matrix because of its relatively low effective refractive index. A low refrac- 
tive index reduces reflection at the sample boundaries and, thereby, reduces the 
size of the diffuse spot and the number of independent scattering channels. 
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In conclusion, for most experiments the requirements on the samples are not very 
stringent and allow for a wide range of materials to be used. We successfully applied 
our wavefront shaping method to daisy petals, porous gallium phosphide, TiOz pig- 
ment, ZnO pigment, white airbrush paint, eggshell, stacked layers of 3M Scotch tape, 
paper and even a baby tooth. 

2.4.1 Sample preparation method 

We developed a spray painting technique to fabricate layers of strongly scattering 
material. The technique allows the fabrication of flat, homogeneous layers with a 
thickness of around 5 pm and more. Moreover, it is very easy to use different materi- 
als or to dope the sample with fluorescent markers. We first give the recipe for ZnO 
samples with embedded fluorescent nanospheres. These samples were used for the 
experiment in Chapter 5. Then we explain how different samples were made using 
the same technique. 

1. Substrate cleaning Standard 40mm x 24mm microscope cover slips with a thick- 

ness of 160 pm were used as a substrate. The substrates were first rinsed with 
acetone to remove any residual organic material. Then they were rinsed with 
isopropanol, a solvent that is known to leave no drying stains, and left to dry. 

2. Paint preparation First, the nanosphere suspension^ was diluted by a factor of 10^ 

(2.5 pi suspension in 250ml water). Then, a suspension was made by mixing 
2.5 g of ZnO powder^ with 7.3 ml water. The suspension was stirred on a roller 
bank for 30 minutes and then placed in an ultrasonic bath for 15 minutes. Fi- 
nally, 0.73 ml of the diluted nanosphere solution was added. The suspension 
was again stirred for 30 minutes and then placed in the ultrasonic bath for 15 
seconds. 

3. Spray painting The paint was sprayed onto the substrate with an airbrush^". The 

airbrush was operated at an air pressure of 2.3 bar. The paint was sprayed from 
approximately 20 cm distance to allow for a homogeneous coverage. The empty 
substrates were taped to an underground with an inclination of approximately 
45°. Spraying covered these substrates with a thin wet film of paint. Direcdy 
after spraying, the samples were left to dry horizontally for two hours. 

The resulting samples are homogeneous, flat up to a variation of 1 pm, and opaque 
(see Fig 2.13). The adhesion between the ZnO and the glass is remarkably good, es- 
pecially since no binding agent was used at all. The thickness of the samples was de- 
termined around a scratch in the sample surface using optical microscopy or Dektak 
profilometry. At A = 532 nm, the mean free path is 0.7±0.2 pm, which was determined 
from total transmission measurements (see Section 5.A.3). 



^ Duke Scientific red fluorescent nanosplieres. Diameter 0.30 ± 5%(ji,m. Suspension witii 1% solids in 
water. 6.7 ■ 10^^ spheres / ml. Dyed with FireFly^'^ excitation maximum 542nm, emission maximum 
612 nm. 

^ Sigma-Aldrich Co. Zinc Oxide powder, < 1 (im, 99.9% ZnO. With a scanning electron microscope (SEM) 

the average grain size was determined to be 200 nm 
^"Evolution Solo Airbrush from Harder & Steenbeck, 0.4 mm needle diameter. 
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Figure 2.13: a) Scanning eiectron microscopy image of a ZnO sample (side view). Botli tine sample-air 
interface (top) and the sample-glass interface (bottom) are visible, b) Pfiotographi of a 
ZnO sample (on a background of millimeter paper). 



Various other suspensions were used to make samples. The most homogeneous 
samples were obtained with a mixture of 2.5 parts airbrush paint^^ to 1 part acrylic 
medium^^. We also created samples using a suspension of 5 g of rutile TiOz pigment^^ 
in 10 ml airbrush thinner. Both the samples with airbrush paint and the samples with 
rutile pigment are weakly fluorescent, which makes it hard to distinguish between 
the probe signal and the background fluorescence. Therefore, for the experiments in 
Chapter 5 we used ZnO samples, which do not fluoresce in the relevant wavelength 
area. 



2.5 Control program 

All elements of the experiment are controlled and synchronized by a control program. 
This program manages frame grabbing, detector triggering, video hardware accelera- 
tion, optimization algorithms, user input, realtime visual feedback, data storage and 
more. To keep a program of this scope maintainable we used a component based 
development strategy. Instead of thinking of the program as a whole, we separated 
it into self-contained components that offer a predefined service. An excellent intro- 
duction into the ideas behind of component based programming and agile software 
development can be found in Ref. (author?) 26. 

Each component implements one or more interfaces. An interface is a specifica- 
tion of a set of properties (publicly accessible data) and methods (executable func- 
tions) with a certain functionality. For example, the interface of a camera component 
defines a method for triggering the camera, and it defines properties such as the shut- 

"Hansa pro-color opaque white. Pigmented, water-based acrylic polymer paint. Pigment: Ti02 PW 6, 

CI#77891 (specifications from manufacturer material safety data sheet). 
'^Schminke 50602 acrylic airbrush thinner 
i^Sachdeben Rutil R210 #1060377/040 
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ter time and the region of interest. Interfaces are defined in the most general terms 
possible and hide all details of the implementation. 

The only way to interact with a component is through its external interface. This 
gives the possibility to create a different component that implements the same inter- 
face with different underlying details, for example using a different low-level camera 
driver. From the outside, these two components are indistinguishable. This means 
that we can exchange these two components transparently. [27] 

The interchangeability of components also makes it possible to extend the func- 
tionality of code without changing it. This tantalizing idea is known as the Open- 
Closed principle [28]. For example, the first algorithm we wrote for shaping a wave- 
front was developed with a photodiode detector and a phase-mostly modulation con- 
figuration. Because of the abstraction that was used, we were able to extend this al- 
gorithm to work with an EMCCD camera and use the phase-amplitude modulation 
scheme simply by connecting two different components to the algorithm. Not a sin- 
gle line in the code of the algorithm was changed. 

2.5.1 LabView and C++ mixed programming 

We chose a mixed programming environment for developing the control program. 
We used C++ to implement reusable components for the detectors, algorithms and 
modulation schemes. Then, we used LabView to construct the measurement proce- 
dure for a specific experiment. The rationale behind this separation to make writing 
a measurement procedure as simple as possible. 

LabView is a visual scripting environment that allows rapid development of simple 
measurement procedures. Moreover, it is very easy to create a user interface for set- 
ting the parameters of the experiment and for giving dynamic feedback in the form of 
graphs and false-color images. However, LabView is not suitable for creating complex 
programs or for performing time critical operations. Also, LabView code usually has 
a limited reusability. 

C++, on the other hand, is a programming language that provides a low level con- 
trol over the hardware and allows for writing time critical code. Moreover, C++ is a 
structured programming language that offers many mechanisms for code reuse and 
abstraction. Additionally, the Microsoft Visual Studio programming environment al- 
lows easy navigation through a complex project and contains various useful tools 
such as an excellent debugger. 

2.5.2 Component Object IVIodel 

The Component Object Model (COM) is a standard for writing reusable components. 
A very accessible introduction to COM is given in Ref (author?) 29. COM was intro- 
duced by Microsoft to enable interprocess communication and dynamic object cre- 
ation in a way that is independent of the programming language that is used. Nowa- 
days, COM and ActiveX (one of the extensions to COM), are supported on various 
operating systems. 

COM/ActiveX components are listed in the Windows registry. The registry contains 
information about the location of the components executable code, as well as a de- 
scription of the interface(s) of the component. In almost any modern programming 
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Figure 2.14: Simplified class diagram for the most important components of the control program. 
Blocl<s indicate classes. The top part of a box contains the class name, the middle part 
contains its properties and the lower part its methods. Arrows denote inheritance. E.g. 
AVT camera is an implementation of IDetector, which, in turn, derives from Rootclass. 



language, including Lab View, it is possible to select a component from the registry 
and use it as if it was written in that programming language. It is even possible cre- 
ate a COM object on a remote computer and access it as if it is a local object. The 
operating system makes sure that all function calls are transparently relayed over the 
network connection. We successfully used this feature to access the Hamamatsu EM- 
CCD camera that was operated by a dedicated computer. 



2.5.3 Global structure of the control program 

We defined four interfaces for describing different functions in the control program. 
These interfaces are 'IDetector' (for accessing detectors and cameras), TModulator' 
(for programming the light modulator), 'lAlgorithm' (that implements a wavefront 
shaping algorithm) and TSample' (for storing and transferring measurement data). 
All interfaces have some common functionality that is embedded in a common root 
class. A simplified class diagram of the interface hierarchy is shown in Fig. 2.14. We 
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first discuss thiese interfaces and tlien explain how these components are connected 
to form a measurement program. 

Rootclass (IDataAccess) We constructed a common base interface for all compo- 
nents. It defines functions for accessing the component's data in a uniform 
way. All publicly accessible data of a component is organized in named 'chan- 
nels'. Access to the data is, in many ways, independent of the data type. For 
example SetChannelValue("Shutter", 1000) sets the shutter time of a camera to 
1000 ^Ls and SetChannelValue("Phase", 0) sets the phase of all segments on a 
modulator to 0. Some data channels represent connections to a different com- 
ponent (more on this below). 

The property Componentlnfo contains a text string with the description of the 
component, including the names and values of all data channels as well as the 
information of all connected components. All measurement programs save 
this string to the disk so that the values of all parameters in the experiment 
are stored. 

■Detector and ISample The detector interface encompasses everything that returns a 
signal. Its channels correspond to settings of the detector, such as the shutter 
time for a camera. When StartMeasurement() is called, the detector returns a 
sample object containing all measured data. Data acquisition happens asyn- 
chronously: the sample object is returned before the measurement and the 
data processing is complete. Only when GetChannelValue() is called on one 
of the data channels of the sample the computer waits until that data is avail- 
able. It is also possible to explicitly wait for the measurement to finish by calling 
Wait() on the sample. 

The detector has a property MeasurementDuration, which is the maximum time 
that the measurement will take. This property is used for synchronizing the 
modulator and the measurements. 

There are many implementations of the IDetector interface. We wrote a com- 
ponent for each of the types of cameras and for the National Instruments data 
acquisition hardware (Nl-daq). We also wrote image processing components 
(ImProc) that analyze a camera image and return a scalar value. Such 'virtual 
detector' components are completely interchangeable with real detectors since 
they have the same interface. 

■Modulator The modulator interface standardizes a light modulator. A modulator 
component has channels that contain matrices for the phase and, depending 
on the modulation scheme, amplitude. The properties SettleTime and IdleTime 
affect the synchronization of the measurements (see Section 2.2.1). The func- 
tion RefresliO is called to display the new image after the phase and/or ampli- 
tude channels have been changed. 

For testing and characterization purposes, different types of modulator com- 
ponents have been created, such as the one that displays a Ronchi grating. 

■Algoritlim The algorithm component contains code for the feedback procedure. The 
three different algorithms that are listed in Fig. 2.14 are discussed in Chapter 7. 
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Figure 2.15: Structure of a typical wavefront optimization program. Rectangles, ActiveX components 
written in C++ are instantiated in a LabView environment and connected to each other; 
lines, connection between a component and a channel (diamond) on a different compo- 
nent; ellipse, a user interface constructed in LabView program provides feedback during 
optimization. 



The algorithm is operated by connecting a detector and a modulator object and 
then repeatedly calling lterate(). 

The LabView program instantiates the required objects and sets parameters such as 
the shutter time, the number of segments, the position of the target, etc. Then, the 
program connects the components together. The structure of a typical wavefront op- 
timization program is shown in Fig. 2.15. A modulator and an image processing node 
are connected to the algorithm. The image processing node takes two inputs, one 
for the camera and one for a reference detector. The node averages the camera im- 
age over the target area and divides this value by the reference signal to compensate 
for variations in the incident intensity. The image processing node also takes care 
that the camera and the reference detector are triggered simultaneously. A LabView 
user interface displays the current status of the algorithm {percentage done, optimal 
wavefront settings measured so far, etc.). 

It proved to be relatively easy to quickly build a program with these predefined 
components. Moreover, most of the LabView code could also be reused in different 
experiments by simply exchanging components. The design is open for extension; 
different cameras, algorithms, target functions and modulation schemes were added 
without affecting other parts of the program. All in all, the flexible design helped in 
quickly setting up an experiment. 



2.6 Conclusions and outlook 

In this chapter, we described the setup that was used for our wavefront shaping ex- 
periments. The complete apparatus was designed from scratch. We discussed the 
considerations that played a role in designing the experiment. Especially the wave- 
front synthesizer turned out to be more involved than the cartoon in Fig. 1.1 suggests. 
We first used a state-of-the-art method to achieve phase-mostly modulation with a 
commercial liquid crystal display. Then, we developed a new method that allows full 
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control over both phase and amplitude of the light. We also discussed a versatile fab- 
rication method that we developed for making samples. This method allowed us to 
make samples with different pigments and to embed fluorescent nanospheres. 

We wrote a control program that executes the optimization algorithm and handles 
all data flow in the experiment. The program was designed using a component based 
development strategy. All detailed programming code was packaged in reusable and 
interchangeable components. Therefore, it was easy to quickly try alternative mea- 
surement schemes, to test new algorithms, and to replace detectors and cameras in a 
late stage of the project. 

The setup that is discussed in this chapter is capable to create an optimal wave- 
front in a matter of minutes. Although this is an impressive achievement by itself, 
our method cannot yet be applied to media with rapidly moving scatterers such as 
perfused skin or turbid fluids. We believe that it is possible to speed up the mea- 
surements considerably. Currently the speed is limited by the response time of the 
LCD (~ 50 ms). However, recently adaptive optics systems for turbulence correction 
have already been demonstrated to reach a correction speed of 800 Hz. [30] Modu- 
lators that use ferroelectric liquid crystals (PLC) or transparent lanthanum-modified 
lead zirconate titanate (PLZT) as a switching medium can have a response time that 
is as short as 3 \is. Even faster switching is achieved with micromachined mirrors; 
such light modulators were shown to have a response frequency of up to 10 MHz. [31] 
Typical frequencies for a centimeter of living human tissue are in the order of 1 — 
100 kHz[32, 33] and we believe that these rates can be achieved with our method in 
the near future. 
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Chapter 3 

Focusing coherent light through 
opaque strongly scattering media 



Abstract: We report focusing of coherent light through opaque scattering materials, 
by control of the incident wavefront. The multiply scattered light forms a focus that 
is up to a factor 1000 brighter than the normal diffuse transmission. 
[This chapter has been published as: I. M. Vellekoop and A. P. Mosk, Opt. Lett. 32, 
2309-2311 (2007)] 



Random scattering of light is what makes materials such as white paint, milk, or hu- 
man tissue, opaque. In these materials, repeated scattering and interference distort 
the incident wavefront so strongly that all spatial coherence is lost.[l] Incident co- 
herent light diffuses through the medium and forms a volume speckle field which 
has no correlations on a distance larger than the wavelength of light. The complete 
scrambling of the field makes it impossible to control light propagation using the well 
established wavefront correction methods of adaptive optics (see e.g. [2]). 

We demonstrate focusing of coherent light through disordered scattering media 
by the construction of wavefronts that invert diffusion of light. Our method relies 
on interference and is universally applicable to scattering objects regardless of their 
constitution and scattering strength. We envision that, with such active control, ran- 
dom scattering wiU become beneficial, rather than detrimental, to imaging[l] and 
communication[3, 4, 5]. 

Figure 3.1 shows the principle of the experiment. Normally, incident light is scat- 
tered by the sample and forms a random speckle pattern (Fig. 3.1a). The goal is to 
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Figure 3.1 : Design of the experiment, a) A plane wave is focused on a disordered medium, a specl<le 
pattern is transmitted, b) Tfie wavefront of the incident light is shaped so that scattering 
makes the light focus at a predefined target. 
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Figure 3.2: Schematic of the apparatus. A 632.8 nm HeNe laser beam is expanded and reflected 
off a Holoeye LC-R 2500 liquid crystal spatial light modulator (SLM). Polarization optics 
select a phase-mostly modulation mode. The SLM is imaged onto the entrance pupil of 
the objective with a 1 :3 demagnifying lens system (not shown). The objective is overfilled; 
we only use segments that fall inside the pupil. The shaped wavefront is focused on the 
strongly scattering sample (S) and a CCD-camera images the transmitted intensity pattern. 
A/4, quarter wave plate. A/2, half wave plate. M, mirror. BS, 50% non-polarizing beam 
splitter. P, polarizer. 



match the incident wavefront to the sample, so that the scattered light is focused in a 
specified target area (Fig. 3.1b). 



3.1 Experiment 

The experimental setup for constructing such wavefronts is shown in Figure 3.2. Light 
from a 632.8 nm HeNe laser is spatially modulated by a liquid crystal phase modulator 
and focused on an opaque, strongly scattering sample. The number of degrees of 
freedom of the modulator is reduced by grouping pixels into a variable number [N] 
of square segments. A charge coupled device (CCD) camera monitors the intensity 
in the target focus and provides feedback for an algorithm that programs the phase 
modulator. 

We performed first tests of inverse wave diffusion using rutile (Ti02) pigment, 
which is one of the most strongly scattering materials known. The sample con- 
sists of an opaque, 10.1 jim thick layer of rutile[6] with a transport mean free path 
of 0.55 ± O.lOtxm measured at A = 632.8 nm. Since in this sample the transmitted 
light is scattered hundreds of times, there is no direct relation between the incident 
wavefront and the transmitted image. [7, 8] 

In Fig. 3.3 we show the intensity pattern of the transmitted light. In Fig. 3.3a we 
see the pattern that was transmitted when a plane wave was focused onto the sam- 
ple. The light formed a typical random speckle pattern with a low intensity. We then 
optimized the wavefront so that the transmitted light focused to a target area with 
the size of a single speckle. The result for a wavefront composed of 3228 individually 
controlled segments is seen in Fig. 3.3b, where a single bright spot stands out clearly 
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Figure 3.3: Transmission tlnrougln a strongly scattering sample consisting of Ti02 pigment, a) Trans- 
mission micrograph with an unshaped incident beam, b) Transmission after optimization 
for focusing at a single target. The scattered light is focused to a spot that is 1000 times 
brighter than the original speckle pattern, c) Multi-beam optimization. The disordered 
medium generates five sharp foci at the defined positions. Figures 3.3a to 3.3c are pre- 
sented on the same logarithmic color scale that is normalized to the average transmission 
before optimization, d) Phase of the incident wavefront used to form Fig. 3.3c. 



against the diffuse background. The focus was over a factor 1000 more intense than 
the non-optimized speckle pattern. By adjusting the target function used as feedback 
it is also possible to optimize multiple foci simultaneously, as is shown in Fig. 3.3c 
where a pattern of five spots was optimized. Each of the spots has an intensity of ap- 
proximately 200 times the original diffuse intensity. In Fig. 3.3d we show the phase of 
the incident wavefront corresponding to Fig. 3.3c. Neighboring segments are uncor- 
related, which indicates that the sample fully scrambles the incident wavefront. 



3.2 Algorithm 

The algorithm that constructs the inverse diffusion wavefront uses the linearity of the 
scattering process. The transmitted field in the target, £,„, is a linear combination of 
the fields coming from the N different segments of the modulator. 



N 



(3.1) 
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Figure 3.4: Measured intensity enhancement as a function of ttie number of segments. Squares, 
sample in focus, triangles, sample 100 |im behind focus, solid line. Ideal enhancement 
(Eq. (3.2)), dotted line, corrected for residual amplitude modulation and finite persistence 
time of T„ = 5400s. The experimental uncertainty is of the order of the symbol size. 



where j4„ and ^„ are, respectively, the amplitude and phase of the light reflected from 
segment n. Scattering in the sample and propagation through the optical system 
is described by the elements ?„„ of the unknown transmission matrix. Clearly, the 
magnitude of E,„ will be the highest when all terms in Eq. (3.1) are in phase. We 
determine the optimal phase for a single segment at a time by cycling its phase from 
to 2n. For each segment we store the phase at which the target intensity is the highest. 
At that point the contribution of the segment is in phase with the already present 
diffuse background. After the measurements have been performed for all segments, 
the phase of the segments is set to their stored values. Now the contributions from all 
segments interfere constructively and the target intensity is at the global maximum. A 
pre-optimization with a small number of segments significantly improves the signal 
to noise ratio. This method is generally applicable to linear systems and does not rely 
on time reversal symmetry or absence of absorption. Although mathematically this 
algorithm is the most efficient, in noisy experimental conditions adaptive learning 
algorithms [9] might be more effective and an investigation in such algorithms is on 
its way (also see Chapter 7). 



3.3 Scaling of the enhancement 

The maximum intensity enhancement that can be reached is related to the num- 
ber of segments that are used to describe the incident wavefront. For a disordered 
medium the constants tmn are statistically independent and obey a circular Gaussian 
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distribution[8, 10, 11, 12] and the expected enhancement r] - defined as the ratio be- 
tween the optimized intensity and the average intensity before optimization - can be 
calculated, 

r]^^iN-l) + l. (3.2) 

It was assumed that all segments of the phase modulator contribute equally to the to- 
tal incident intensity We expect the linear scaling behavior to be universal as Eq. (3.2) 
contains no parameters. Also, since we are free to choose the basis for Eq. (3.1), we 
expect to find the same enhancement regardless of whether the target is a focus or 
a far-field beam and regardless of how the shaped wavefront is projected onto the 
sample. Interesting correlations between the transmission matrix elements, which 
will cause corrections on Eq. (3.2), are predicted when N approaches the total num- 
ber of mesoscopic channels. [11, 12] With our current apparatus we are far from this 
regime and no deviation from Eq. (3.2) is expected. 

We tested the universal scaling behavior implied by Eq. (3.2) by changing N. Using 
the same Ti02 sample as before, the algorithm was targeted to construct a coUimated 
beam. In Fig. 3.4 the enhancement is plotted as a function of the number of segments 
for different focusing conditions. The linear relation between the enhancement and 
the number of segments is evident until the enhancement saturates at 17 = 1000. All 
measured enhancements were slightly below the theoretical maximum. This is un- 
derstandable since all perturbations move the system away from the global maxi- 
mum. The main reason for deviations from the optimal wavefront is residual ampli- 
tude modulation in the phase modulator, which introduced an uncontrolled bias in 
the field amounting to 14% of the total intensity. 

The saturation of the enhancement is the result of slow changes in the speckle pat- 
tern. This instability effectively limited the number of segments for which the optimal 
phase could be measured. We estimate that the effective enhancement decreases to 
i7eft — ri/il + NT/Tp), where T = 1.2s is the time needed for one measurement and 
the persistence time Tp — 5400s is the timescale at which the speckle pattern of the 
Ti02 sample remains stable. Depending on the environmental conditions, Tp can be 
considerably higher and enhancements of over two thousand have been measured 
overnight. 

To verify the universal applicability of inversion of wave diffusion, we used a variety 
of materials of natural origin. Table 3.1 lists the intensity enhancement for different 
materials we used. Although the samples vary in thickness, composition and scat- 
tering strength, they were all able to focus a properly prepared wavefront to a sharp 
spot. The intensity enhancement varies between 60 and 1000. The main reason for 
this variation is that the persistence time is not the same for all materials. 



3.4 Conclusion 

In summary, our results show that precise control of diffuse light is possible using 
an optimal, non-iterative algorithm; light can be directed through opaque objects to 
form one or multiple foci. The brightness of the focal spot is explained by a model 
based on statistical optics. We expect inverse wave diffusion to have applications 
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Sample L{\i.m) rj{±15%) N 

RutileTiOz 10.1±0.3 1080 3228 

Daisy petal, fresh 43±5 64 208 

Daisy petal, dried 37±5 630 1664 

Chicken eggshell 430±30 250 3228 

Human tooth 1500±100 70 208 

Table 3.1 : Measured intensity enhancement for different materials. L, sample thickness (± surface 
roughness); rj, maximum enhancement reached; Af, number of segments used by the algo- 
rithm to describe the wavefront. 



in imaging and light delivery in scattering media, possibly including metal nano- 
structures[13]. Dynamic measurements in biological tissue are possible when the 
time required for achieving a focus can be reduced to below 1 ms per segment[14, 
15]; we estimate that this timescale is technologically possible with the use of fast 
phase modulators[16]. Furthermore, the high degree of control over the scattered 
light should permit experimental verification of random matrix theories for the trans- 
port of light[ll, 12]. 
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Chapter 4 

The focusing resolution of 
opaque lenses 



In this chapter, we investigate both experimentally and theoretically the focusing res- 
olution of the opaque lenses that were introduced in Section 1.1. For every optical 
system there is a limit on how sharp it can focus light. The smallest focus that can the- 
oretically be achieved with an aberration free lens is called the diffraction limit. The 
diffraction limit is determined by the aperture and focal distance of the lens. Here 
we examine how these concepts translate to the situation of an opaque lens. We use 
the experimental geometry that is shown in Fig. 4.1 to measure the focal v«dth of an 
opaque lens and then compare this width to the theoretical diffraction limit of an 
ordinary lens. 

In Section 4.1 we explain how scattering inside an opaque lens helps focusing. The 
focusing resolution of an opaque lens is determined experimentally in Section 4.2. In 
Section 4.3 we compare the shape of the focus with the measured speckle correlation 
function. To theoretically understand the focusing capacity of opaque lenses, a ma- 
trix model does not suffice. Therefore, we developed a continuous field description 
of opaque lenses in Section 4.4. This new theoretical framework explains the experi- 
mental results. Moreover, the theory predicts the existence of the optical equivalent 




fi = 200 mm 



Figure 4.1: Schematic of tiie experimental geometry. Liglit from a phiase modulator is imaged on the 
center plane of a lens (imaging telescope not shown). The numerical aperture of the lens is 
controlled by a pinhole. A CCD camera is positioned in the focal plane of the lens. Between 
the camera and the lens a strongly scattering sample is placed. The sample can be moved 
to change the distance to the camera. 
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Figure 4.2: Scattering helps wavefront shaping and focusing. A single incident beam contributes to all 
transmitted wave vectors. As a result, a sharp focus can be generated even if the wavefront 
synthesizer only controls light in a small solid angle. 



of the acoustical time-reversal lenses that focus sound beyond the diffraction limit[l] . 



4.1 Wavefront shaping with an opaque lens 

In practice, it is extremely difficult to design (systems of) lenses that have almost no 
aberrations. [2, 3] A good microscope objective, for instance, contains up to twenty 
lens elements for compensating various aberrations. [4] An alternative approach to 
achieve a diffraction limited focus is to use adaptive optics to compensate for the 
aberrations in the lens system. [5] However, the adaptive optics system that gener- 
ates the incident wavefront may have many limitations, such as a maximum num- 
ber of control segments, residual phase-amplitude cross-modulation, a limited max- 
imal wavefront curvature, etcetera. AH such limitations cause systematical deviations 
from the diffraction limited focus. 

We found that, surprisingly, the shape of the focus of an opaque lens is not affected 
by experimental limitations of the wavefront modulator. Paradoxically, the scattering 
medium makes that we have apparently more control over the scattered field than we 
have over the incident field. A similar effect has been observed in time-reversal ex- 
periments using ultrasound and microwaves: a limited number of broadband trans- 
mitters suffices to generate a diffraction limited focus through a disordered medium 
(see [6] for an excellent review on time-reversal focusing). Although our experiments 
are performed at a single (optical) frequency, many of the concepts developed for 
time-reversal experiments are still valid. 

In Fig. 4.2, we illustrate how scattering improves focusing. Disordered scattering in 
the sample randomly redistributes the incident light over all scattering channels (see 
e.g. [7, 8, 9]). In a microwave experiment, it was demonstrated that even evanescent 
channels are excited.[10] Scattering makes that light from a single segment of the in- 
cident beam contributes to all wave vectors of the desired field. By tuning the relative 
phases of all segments, their contributions to the desired field add up constructively. 
The more segments are used, the better the desired field can be approximated. The 
constructed field will never be perfect as the number of degrees of freedom in the 
transmitted field cannot exceed the number of control segments. However, thanks 
to the disordered scattering any deviation in the generated field is, on average, uni- 
formly distributed over all outgoing channels, resulting in a low overall intensity con- 
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tribution. Imperfections in the generated wavefront cause a loss of contrast with the 
speckle background. 

4.2 Measured focusing resolution of an opaque lens 

The experimental geometry used for analyzing the focus of an opaque lens is shown 
in Fig. 4.1. The sample is a 6p,m layer of airbrush paint on a standard microscope 
cover slip. The wavefront illuminating the sample is shaped to increase the intensity 
on a single pixel of the camera. The algorithm that shapes the wavefront is the same 
as was used in Chapter 3; it does not take into account the shape of the focus. The 
distance between the sample and the camera can be varied to analyze the focusing 
resolution of the opaque lens at different focal lengths. 

We compare the measured focal width to the calculated diffraction limit. The op- 
timal focus of a homogeneously illuminated lens with a circular aperture is an Airy 
disk.[2] The half width at half maximum of the Airy disk is w — 0.51A/NA, where NA 
is the numerical aperture of the lens. This radius is almost equal to the Abbe diffrac- 
tion limit[ll] of If = A/(2NA). To calculate the diffraction limit of the opaque lens, we 
need to determine its numerical aperture NA2 = 02/(2/2), where D2 is the diameter 
of the illuminated area on sample and /2 is the focal length of the opaque lens. 

The diameter of the spot on the sample depends on what wavefront is generated. 
When the glass lens in front of the aperture is illuminated with a plane wave, the di- 
ameter of the spot on the sample equals D2 — 2/2 NAi, where the numerical aperture 
of the glass lens is NAi = Di/(2/i) — 5 • 10 ^ with Di — 2.1 mm the diameter of the 
aperture and /i — 200 mm the focal length of the lens. 

When the modulator is programmed to generate a shaped wavefront, diffraction at 
the lens aperture will cause the beam to diverge. Segments of the phase modulator 
are imaged onto the aperture to squares with sides of length Ds — 0.03 mm. The angle 
of diffraction for a square source equals a = {4X)/{7i^Ds) = 0.008 rad (see e.g. [12]). 
Combining the effects of the focusing and diffraction, we find the spot diameter 

D2 = 2/2NAi+2(/i-/2)a, (4.1) 

where / — 200 mm is the focal length of the lens. The radius (half width at half maxi- 
mum intensity) of the diffraction limited focus at the detector plane is given by 

Af2 

■" (4.2) 



2/2NAi + (/i-/2)a 



The experimental results are displayed in Fig. 4.3. All widths are relative to the 
diffraction limit of the glass lens {wo = Sljim). Without sample, the sharpest focus 
that was obtained by optimizing the wavefront has a width of 1.2 times this diffrac- 
tion limit. The fact that the diffraction limit was not reached clearly demonstrates 
the limitations of the wavefront synthesizer. With a sample, however, these limita- 
tions do not degrade the sharpness of the focus. The width of the focus is exactly the 
diffraction limit that is predicted by Eq. (4.2). 
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Figure 4.3: Relative focal width w as function of sample-camera distance fi. Diamonds, measured val- 
ues with a sample present; solid curve, theoretical diffraction limit; circle (at fi =200 mm), 
width of the focus after optimizing without sample; dashed line, size of a single pixel of the 
camera, detection limit. All widths are relative to the diffraction limit of the glass lens (dotted 
line, Wa = 31|jm). 



4.3 Measured relation between the focus and the speckle 
correlation function 

It appears that the focus of an opaque lens always has the same size as a typical 
speckle. In this section, we analyze this effect by comparing the measured profile 
of the focus to the measured speckle correlation function. The speckle correlation 
function C is a measure for the shape of a typical speckle. [7] It is defined as 



C(xi,yi;x2,y2); 



(7(xi,yi)J(x2,y2)) 
(7(xi,yi)) (7(x2,y2)) 



1, 



(4.3) 



where the brackets denote spatial averaging over all speckles. We experimentally in- 
vestigated the relation between the speckle correlation function and the generated 
focus. First, we used the geometry as shown in Fig. 4.1 with the sample as close to the 
camera as possible (26 mm). The resulting spot has a radius of 6.6 ± 3 p,m, which is 
hardly larger than the size of a single pixel on the camera (6.45 [xm) . The speckle cor- 
relation function and the intensity profile of the spot overlap exactly, as can be seen 
in Fig. 4.4. 

The same experiment was performed in a different geometry and with different 
samples. This time the incident light was focused using a microscope objective and 
the transmitted light was imaged using a second objective (see Chapter 3) . For a rela- 
tively thick sample (egg shell, 430 ± 30 ^xm thick) we placed the detecting objective at 
a distance of 3.5 ± 0.5 mm and found a very narrow focus of I . I ± 0.2 jim (see Fig. 4.5). 
In Fig. 4.6, the results for a thinner sample (I0.1±0.3 jim of rutile Ti02 pigment) in the 
same configuration is shown. For this sample, we find a vnder focus of 11 ± 0.2 ^xm. In 
both cases, the profile of the focus overlaps with the speckle correlation function. 
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Figure 4.4: Intensity profile of the focus at y = (solid curve) and speckle correlation function (dashed 
curve) for an 8|im thick layer of airbrush paint creating a focus at 26mm distance. The 
speckle correlation function was measured with a randomly generated incident wavefront. 
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Figure 4.5: Intensity profile of the focus at y = (solid curve) and speckle correlation function (dashed 
curve) for an egg shell of 430±30|im thick creating a focus at 3.5 ±0.5 mm distance. 
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Figure 4.6: Intensity profile of the focus at y = (solid curve) and speckle correlation function (dashed 
curve) for a 10.1 ± 0.3 |im thick sample of white paint creating a focus at 3.5± 0.5mm 
distance. The central point of the correlation function was removed since it was affected by 
measurement noise in the camera image. 
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4.4 Continuous field theory for opaque lenses 

After optimizing to create a focus, a fraction of the field at the back of the sample 
forms a converging beam, while the rest of the field forms a random speckle. In Chap- 
ter 3, we introduced a matrix model to explain the intensity of the focus of an opaque 
lens. Although that model proved suitable for describing the intensity in the focus of 
an opaque lens, it cannot predict the shape of the focus. We now develop a continu- 
ous field formalism that explains both the intensity and the shape of the focus. This 
formalism also explains the similarity of the focus profile and the speckle correlation 
function that was observed in the experiment. 

4.4.1 Continuous field formalism 

Light propagates from the modulator to the sample, then it is transmitted by the sam- 
ple, and finally it reaches the focal plane where the target detector is placed. We sep- 
arate light propagation into two parts. The first part starts at the modulator and ends 
at the back of the sample. The second part is the propagation from the back of the 
sample to the focal plane. 

We define coordinates in the focal plane as tb [Zb — fi as defined in Fig. 4.1), on the 
back surface of the sample as vt [Zk — Q], and in the plane of the phase modulator as 
Va . The target for optimization is at point vp . The field at the back of the sample is 
given by 

JV 
E{r!c)=^ta{Tk)Ea, (4.4) 

a 

where the summation is over all segments of the phase modulator. We assume the 
field to be constant within a segment, so that we can write EijCa) — Ea ■ Just behind the 
modulator, light propagation is nearly perpendicular to the surface of the modulator. 
We can, therefore, use the Fresnel-Kirchhoff diffraction formula (see e.g. [13, 14]) to 
find coefficient tairt] 

tairk]=j\\ d2r,G(rfc,rJ. (4.5) 

The integration is over the surface of a single segment of the light modulator (Sa) and 
G{rt,ra) is the unknown Green function for propagating from the modulator to the 
back of the sample. Light propagation from the back of the sample towards the focal 
plane is described by the free space Green function g{rb — rt] 

£(rfc)= d'rkg{rb-rk]jE{rk). (4.6) 

Field propagation from a single segment a to the focal plane is described by a trans- 
mission coefficient ta(rb) 

fjrb)= dhkg{rb-rk]jtairk). (4.7) 

The role of ta{rb] is similar to that of tba in the matrix formalism. However, ta{rb) is a 
continuous function of the spatial coordinate and also gives the shape of the focus of 
the opaque lens. 
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4.4.2 Optimized field 

We consider the case where the incident field is shaped with a modulator that per- 
fectly modulates phase and amplitude. In Chapter 8, we will see that the results for 
imperfect modulation or phase only modulation are equal up to a constant prefactor. 
After optimizing for focusing at point rp, the field generated by the modulator is (see 
also Chapter 8) 

Ea^Eotlirp), (4.8) 

where the tilde sign indicates the value of Ea after optimization. Eo is a constant that 
fixes the average intensity at the modulator to /in, it is defined as 

(4.9) 



JV 

7o(r^) = 7in^|rfl(r^)p. (4.10) 

a 

A physical interpretation of loirp) is that it is the intensity in the focal plane averaged 
over all possible configurations of the light modulator that have an average intensity 
of/in- 

The field at the back of the sample after optimization is calculated by substituting 
Eq. (4.7) into Eq. (4.8) and then substituting the result in Eq. (4.4) 



JV 

£(rfc) = y f„(rt) I £o 1 1 d'r[girp-r',)^taW,) 1 . (4.11) 



:^ar,) £o dXg(r^-i4)^, 



The term fa(ri:) t*(r'j.) fluctuates rapidly with the distance Ar^r'j. —rt- The ensemble 
average of this product is the field-field correlation function[15] 

{taivkKir'^)) = (|f«(rfc)|2) ^^^^^^exp(-|Ar|/(2£)), (4.12) 

where fco = 2n/?. and |ta(rfc)|2 is the intensity transmission coefficient from segment 
a to point rj,. Its average value {\ta{rkW) depends on the geometry of the experiment 
and the thickness of the sample. Because the correlation function Eq. (4.12) is very 
sharp, we can take g*(r^ ~'"'/t^ ^"^^ of the integral in Eq. (4.11). Integrating Eq. (4.12) 
over r'j. gives. 

jj dh[{UrkKir',]) = (|f,(r,)p) 2n^^^^^^ ^ (|f.(r,)|2) g, (4.13) 

where we used that 4{kolY "> 1. This assumption is valid even for very strongly scat- 
tering samples. The area given by X^l(2n) corresponds to the size of a single scatter- 
ing channel (see e.g. [16]). The average value of the optimized field at the back of the 
sample now reduces to 

(£(rfc)) = g*{Tfi - Vk)Coh{Tk), (4.14) 
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Figure 4.7: Better than diffraction limited focusing with an opaque lens, a) An ordinary lens at best 
focuses light to the diffraction limit given by the diameter of the incident beam, b) An 
opaque lens focuses light more sharply than the diffraction limit because diffusion in the 
lens increases the diameter of the transmitted converging beam, c) Opaque wafer lens 
made out of gallium phosphide. Reflection at the GaP-air interface fold back the optical 
paths and improve the focusing resolution even further. 



where Co = iXEol{2Tzh^) and /o(r*:)= /mXa \ta{rk)V, just as in Eq. (4.10). 

In the geometry as shown in Fig. 4.1, Ia(tk] is constant in a disk with a diameter D 
given by Eq. (4.1) and zero outside of the disk. In this case, Eq. (4.14) becomes 



[E{rk)) = g*[rp - Vk)CoIoCViCD(rk), 



with 



circo(rj:)E 



fo x2+y2>DV4 
I 1 otherwise. 



(4.15) 



(4.16) 



Equation (4.15) tells us that, up to a constant proportionality factor, the average field 
at the back of the sample is the complex conjugate of the free space propagator. In 
other words, the transmitted field is the same as if the sample time-reversed a wave 
coming from the target focus. Hence, the transmitted light focuses on the target as 
well as is physically possible for a system with an aperture given by chcoivk)- 

In a different geometry, for example when the incident light is focused sharply on 
a thick sample, the intensity distribution /o(ryt) will be a smooth function of the po- 
sition. In such a geometry, an opaque lens differs from an ordinary lens. With an 
ordinary lens the amplitude profile at the back of the lens is equal to the amplitude 
profile of the incident illumination. For an opaque lens, however, the amplitude of 
the controlled field at the back of the sample is proportional to the intensity distribu- 
tion at the back of the sample. In the following section we will see how to make use of 
this difference. 



4.4.3 Better than diffraction limited focusing 

In time-reversal microwave experiments it was observed that scattering improves the 
focusing resolution beyond the diffraction limit by 'folding' the propagation paths in 
a waveguide. [1] In this section, we show that a similar dramatic improvement of the 
focusing resolution can also be obtained in an optical system with a slab geometry. 

With an ordinary lens, the intensity profile of the transmitted is the same as that of 
the incident light. For a given diameter of the incident beam, the lens at best focuses 
to the diffraction limit A//D (see Fig. 4.7a). 
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With an opaque lens, however, diffusion in the lens makes that the transmitted spot 
is larger than the incident beam (Fig. 4.7b). Therefore, an opaque lens focuses light 
more sharply than the best possible ordinary lens. The intensity distribution loirt) at 
the back of the sample depends on the thickness and mean free path of the sample, as 
well as on the extrapolation ratios that describe the boundary conditions. The thicker 
the scattering medium, the sharper the focus. 

The benefit of a thick scattering medium was observed in the experiment. Under 
the same experimental conditions, the focus of a 430 [xm thick opaque lens (Fig. 4.5) 
was ten times narrower than the focus of a 10.1 ^im thick opaque lens (Fig. 4.6). 

A proposed special example of an opaque lens is shown in Fig. 4.7c. A thin layer on 
one side of a gallium phosphide (GaP, n=3.3) wafer is made porous. The porous gal- 
lium phosphide scatters the light into a 2n solid angle. However, due to the high index 
contrast at the back of the wafer, most of the transmitted light is internally reflected 
back to the scattering layer, where it is scattered by the porous GaP for a second time. 
By folding the optical paths this way, the spot at the back of the sample increases 
dramatically, and the focusing resolution improves accordingly. 

We now estimate the factor by which a GaP wafer lens improves the diffraction 
limit. The critical angle for total internal reflection at the GaP-air interface is 9c — 
17.6° . For a wafer with a thickness of L, the incident spot needs to have a diameter of 
at least 2 tan Be L for the reflected internally reflected light to overlap the original spot. 
We neglect the contribution of light that is reflected at an angle of over 45° and we 
take into account only the first reflection at the GaP-air interface. Even with this pes- 
simistic estimate, a 0.5 mm thick opaque a GaP wafer lens already expands a 0.32 mm 
beam to 4 mm before focusing it. These realistic numbers indicate that a 0.5 mm thick 
wafer can replace a 6x beam expander and a perfectly corrected microscope objec- 
tive. In a way, this wafer is a miniature Cassegrain objective. 



4.4.4 Intensity profile of the focus 

We use the general result in Eq. (4.14) to calculate the intensity distribution in the 
focal plane. Using Eqs. (4.6) and (4.14) we find 

£(rfa)=^ d^Vkg[vb-Vk)g*[vp-Vk)CMvk). (4.17) 

In the paraxial limit (see e.g. [14]), we can approximate 

g(rfc-ri:)g*(r/;-rfc)!%i— j-expl — [xk[xp-Xb) + yk[yp-yb)\ \, (4.18) 

with /2 the focal distance of the opaque lens. We isolated phase factor <1> = jCg — x? + 



^P ^b 



Vp-Vb- 



[E{tb)) = -^Co d2rfc7o(rj;)exp ( — ^ [xk{xp - Xb) + yk{yp - Yb)] J ■ (4.19) 
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In Eq. (4.19) we recognize a 2-dimensional Fourier transform. We define the Fourier 
transformed intensity L/q as 

L/o(ri,)= d^rfc7o(rfc)exp — ^ [xkixp - Xb) + ykiyp - yt)] j ■ (4.20) 

To conform to literature (see e.g. [7]), we adopt the notation 

J2-Uo(ri:)} = t/o(rb). (4.21) 

Now, the intensity in the focal plane is given by expanding Q 

nrb]^—^j^,\^{Idrk)}f, (4.22) 

which confirms the observation made in Section 4.4.2 that the opaque lens perfectly 
focuses a field that is proportional to the intensity distribution at the back of the sam- 
ple. The intensity exactly at the center of the focus equals 

N 
271/1 j loJ^pY 



^("•^)=ld7^l 77:^^' (4.23) 



where Ptot = jj (i^^khiXk) is the average total power that is transmitted through the 
sample. Before optimization, the power is distributed over 27r solid angle. Therefore, 

/o(r/}) = Ptot/(27r/|)and 

I[rp)^NUTp), (4.24) 

which means that the intensity in the target has increased by exactly a factor N. For 
phase only modulation, the enhancement reduces to 1 + (AT — l)7r/4 (see Chapter 8). 
The matrix model used in Chapter 3 predicts the same intensity enhancement. How- 
ever, the new formalism clearly shows that the enhancement needs to be measured 
exactly at the maximum of the focus and not, for instance, averaged over the whole 
speckle. This experimentally important detail was not clear from the matrix model. 
Moreover, the new formalism predicts the exact shape of the focus. 

4.4.5 Connection with speckle correlation function 

When the sample is thicker than a transport mean free path, the field at the back 
surface of the sample is uncorrelated on length scales larger than the wavelength of 
light (also see Eq. (4.12)) and, equivalently, the transmitted light has no preferential 
direction. For such a sample, in the paraxial limit the correlation function (Eq. (4.3)) 
only depends on the coordinate difference Arj, = r^ — r^'. The van Cittert-Zernike 
theorem tells us that, under these conditions, the correlation function is given by[7, 
17, 18] 

C(Axb,Ayfe)=-^|J2-{(7(r,)>}^ (4.25) 

tot 

where /(rfc) is the diffuse intensity distribution at the back surface of the sample and 
^{{I{rk))] is defined as in Eq. (4.21). When the phase modulator is configured to 



Conclusion 63 



generate a random field, (/(r/t)} — (loirk))- In this case, it follows from Eq. (4.22) that, 

7(x,y,/2) = C(x,y)iV -^7o(r^). (4.26) 

Equation 4.26 predicts that the intensity profile of the focus is exactly equal to the 
speckle correlation function (up to a constant prefactor). This identity is striking be- 
cause the intensity in the focus is the square of an amplitude, whereas the correlation 
function is an product of intensities. Our continuous field model explains the results 
of Section 4.3, where this identity was observed experimentally in different geome- 
tries. 

4.5 Conclusion 

We observed experimentally that opaque lenses focus light to a diffraction limited 
spot. This experimental observation is explained by a continuous field model that 
was developed in Section 4.4. The model provides a simple method to predict the 
focusing resolution of an opaque lens for different geometries. For instance, a thin 
opaque lens that is illuminated with a flat-top profile focuses light to a diffraction 
limited spot. Theoretically, the focusing resolution can even be improved beyond the 
diffraction limit by making use of the diffusion in a thick opaque object. Extremely 
sharp focusing can be achieved with a GaP wafer lens where the optical paths are 
folded back by total internal reflection. 

The shape of the focus is exactly the same as the correlation function of the speckle 
that is generated by the opaque lens. This experimentally important observation was 
explained with the continuous field model for opaque lenses. 

An opaque lens focuses light to a diffraction limited spot or even sharper. Addi- 
tionally, opaque lenses are extremely thin and extremely cheap compared to a mi- 
croscope objective. Therefore, disposable opaque lenses might replace conventional 
microscope objectives, for example in applications that require tight focusing in a 
polluting or corroding environment, or when space is tight. 
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Chapter 5 

Demixing light paths inside 
disordered metamaterials 



Abstract: We experimentally demonstrate the first method to focus light inside dis- 
ordered photonic metamaterials. In such materials, scattering prevents light fi-om 
forming a geometric focus. Instead of geometric optics, we used multi-path in- 
terference to make the scattering process itself concentrate light on a fluorescent 
nanoscale probe at the target position. Our method uses the fact that the disorder in a 
solid material is fixed in time. Therefore, even disordered light scattering is determin- 
istic. Measurements of the probe's fluorescence provided the information needed to 
construct a specific linear combination of hundreds of incident waves, which inter- 
fere constructively at the probe. 

[This chapter has been published as: I. M. Vellekoop, E. G. van Putten, A. Lagendijk, 
and A. P. Mosk, Opt. Express 16, 67-80 (2008)] 

Photonic metamaterials are materials that gain their optical properties from their 
nanostructure, rather than from their chemical composition. Their sub-wavelength 
elements interact very strongly with the electromagnetic field and provide unprece- 
dented control over the propagation of light[l, 2]. Exciting perspectives of metama- 
terials include field enhancement in photonic crystal cavities [3], negative refractive 
index optics[4, 5, 6], and sub -wavelength imaging[7, 8]. Perhaps even more amaz- 
ing is the use of metamaterials to cloak objects, making them completely invisible 
for a certain frequency of light. [1, 2] While these perspectives require near-perfecdy 
manufactured metamaterials, intentionally disordered metamaterials also find appli- 
cations such as sub -wavelength focusing[9] and enhanced non-linear conversion[10, 
11]. Furthermore, disordered metamaterials are of fundamental interest because of 
their analogy to electrical conductors [12, 13], their quantum optical properties[14], 
and mesoscopic transport phenomena[15, 16, 17]. Recently it was shown experimen- 
tally that localization induced by disorder can confine and guide light. [18] 

To direct light to a desired position in a metamaterial requires thorough knowledge 
of the propagation of light in these structures. Since metamaterials are strongly pho- 
tonic, even weak manufacturing imperfections strongly influence propagation [19, 
20, 21], which makes it even harder to get light to the right place. For disordered 
metamaterials it is, even in theory, impossible to have a priori knowledge of light 
propagation. To treat the metamaterial as a homogeneous medium and pretend to 
focus somewhere inside (Fig. 5.1a) is an erroneous approach; in reality, the structure 
of the material prevents light from forming any focus at all. 

We demonstrate experimentally that light can be delivered to a target inside a dis- 
ordered photonic metamaterial, or in fact any strongly scattering solid, without hav- 
ing any a priori knowledge of the optical properties of the system. Our approach is 
based on the fact that the structure of the metamaterial mixes light entering through 
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Figure 5.1 : Principle of cliannel demixing. a) Conventional way of illuminating a metamaterial sample 
or device. In order to get ligfit to the area inside the metamaterial indicated by the cir- 
cle, one might try to focus light using an ordinary lens. However, the randomly oriented 
nanoparticles in the material mix the incident scattering channels in a complicated way so 
that light does not travel directly to the target. Incident rays will not converge to a geometri- 
cal focus, b) Our new way of illuminating a specific point inside a metamaterial. The wave 
nature of light is used as an advantage. When the phase delay for each of the incident 
scattering channels is set correctly, the channels demix inside the sample; multi-channel 
interference makes the light focus at the desired point. 



different scattering channels (a scattering channel is an angular or spatial mode of 
the optical field[12]). The mixing process is complex and stochastic due to disorder. 
Even so, reciprocity dictates that there always exists a linear combination of incident 
channels that demix inside the metamaterial to focus at the target. In our experiment 
we find this unique configuration and focus light on a target that is deeply embedded 
in a disordered photonic structure (Fig. 5.1b). 

It is well known that the inverse wave problem is ill posed. In other words, even in 
principle optical measurements cannot be used to reconstruct the three-dimensional 
structure of a photonic metamaterial. Likewise, optical reflection and transmission 
measurements cannot reveal the combination of channels that demix inside the ma- 
terial. 

There exist different methods for focusing waves in the presence of scattering. In 
transparent media where geometrical optics applies, delivery of light can be achieved 
very effectively using a variety of adaptive optics wavefront correction techniques 
(see e.g. Ref. (author?) 22). For example, adaptive optics is commonly used to re- 
duce astronomical seeing[23], to improve the resolution of ophthalmic imaging[24], 
or for femtosecond pulse shaping[25]. Unfortunately, these adaptive optics tech- 
niques break down in nanostructured media where geometrical optics is not appli- 
cable. 

For ultrasound and microwaves time-reversal has proven highly effective [9, 26] in 
refocusing waves. Light waves can be refocused back to their source by phase con- 
jugation[27], or they can be focused behind a scattering region by using feedback 
from an external detector[28]. All these methods require physical access to the target 
position, to place a transceiver, an extensive non-linear medium, or a detector. Inside 
a photonic metamaterial, v«th structures smaller than the wavelength of light, these 
methods cannot be used. 

We have developed an approach that does allow focusing of light inside a meta- 
material. Our method combines an advanced feedback and control system with a 
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Figure 5.2: Experimental setup and samples, a) Simplified schematic of the experiment. A 532 nm 
laser is expanded and illuminates a spatial light modulator (SLM) that spatially modulates 
the phase of the reflected light. The SLM is imaged onto the bacl< aperture of a 63x mi- 
croscope objective that focuses the modulated light on a sample. Fluorescence light is 
collected by the same objective and imaged with an EMCCD camera. A computer drives 
the SLM and analyzes the EMCCD images. PBS, polarizing beam splitter; DM, dichroic 
mirror; F, fluorescence filter. Lenses were omitted from the schematic, b) Geometry of the 
samples. A low concentration of 300-nm diameter fluorescent spheres is dispersed in a 
zinc oxide pigment, deposited on a glass substrate, c) SEM image of a sample. 
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Figure 5.3: Experimental demonstration of channel demixing. a) Fluorescence image of a 300 nm 
sphere, embedded in ZnO pigment at a depth of 9.7 |jm. Conventional illumination after 
a scan to find the position where the fluorescence is maximal; approximately 3 times the 
speckle average intensity. Inset, phase of the incident wavefront (plane wave), b) Our 
channel demixing method. Same sphere and sample position as in a, the phases of the 
incident channels were set to the measured optimal values. The resulting fluorescence 
intensity has increased by a factor of 22 with respect to the average intensity Inset, phase 
of the incident wavefront. Inside the sample this seemingly random wavefront transforms 
into a focus. Intensities are in counts per second. 
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nanoscale fluorescent probe that is situated at the desired location. The power of the 
emitted light is a measure for the intensity at the position of the probe. In our exper- 
iment (Fig. 5.2a), we used a liquid crystal spatial light modulator (SLM) that spatially 
modulates the phase of light coming from a green (532 nm) continuous wave laser. 
The pixels of the SLM are grouped into 640 square segments that are imaged onto 
the back aperture of a microscope objective. In this configuration, at the focal plane 
of the objective each segment corresponds to a solid angle encompassing approxi- 
mately one mesoscopic scattering channel. For each of these incoming channels the 
phase can be set individually using the SLM. The microscope objective illuminates a 
sample that is mounted on a xyz- translation stage. The fluorescence light emitted by 
the probe inside the sample is collected by the same microscope objective and im- 
aged with an electron multiplying charge coupled device (EMCCD) camera placed in 
a confocal reflection configuration. Our samples consist of a white layer of zinc oxide 
(ZnO) pigment on a glass substrate (Fig. 5.2b). The pigment grains have an average di- 
ameter of 200 nm, resulting in a sub -wavelength disordered structure (Fig. 5.2c). The 
pigment layer is up to 32 [im thick, and all samples are fully opaque. The transport 
mean free path was measured by fully analyzing transport of light through a series 
of samples (see Section 5.A). The samples were found to be strongly photonic with a 
mean free path of 0.7 ± 0.2 pm, comparable to the wavelength of light. Fluorescent 
spheres with a diameter of 300 ± 15 nm are sparsely dispersed through the sample. 

Conventional focusing of light does not work in disordered photonic media. The 
absence of a ballistic focus was confirmed by performing an extensive 3-dimensional 
scan of the sample with respect to the microscope objective (see Section 5.A). The 
only change in the fluorescence emission was the result of the volume speckle field in 
the medium that is formed by scattered monochromatic light. As the sample moves, 
the weak fluorescence emission fluctuates randomly with the sample position. We 
consequently find that the maximal emission found by scanning the sample is only 
approximately 3 times the typical emission. In numerical simulations we find a sim- 
ilarly weak enhancement (Section 5.B). The fluorescence emission with the sample 
positioned at the brightest speckle is shown in Fig. 5.3a. 

We now cycle the relative phase of a single incident channel between and 2?! 
while keeping the phase of all other channels constant. The intensity at the probe 
is the result of interference between the modulated channel and the other channels, 
and therefore has a sinusoidal dependence on the phase. From the phase and the 
amplitude of the sinusoid we obtain the complex propagation amplitude from the 
modulated channel to the point inside the metamaterial. The complete set of prop- 
agation amplitudes contains all information that is needed to calculate the propaga- 
tion of any incident wave towards the target. The optimal configuration of incident 
channels is found by taking the complex conjugate of the propagation amplitudes, so 
as to fully compensate for the mixing of the channels that will occur. 

Figure 5.3b shows the strong increase in fluorescence emission that was achieved 
with channel demixing. By constructing the optimal configuration of incident chan- 
nels, the intensity distribution inside the sample changed from random speckle to a 
pattern specially tailored for focusing light at the target. We find that the amplitudes 
for neighboring channels were uncorrelated, which is an indication that the channels 
were originally fully mixed. Light from different incident channels mixed while prop- 
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Figure 5.4: Measured intensity for targets at different depths inside a 29±3|im thicl< sample. Stars, 
results of our channel demixing technique, starting from a random sample position. An 
intensity was reached that is on average 20.4 times higher than the diffuse background 
(dotted line). Circles, results of our channel demixing technique starting on a bright speckle, 
demonstrating that starting from a bright spot is not advantageous. Squares, maximum 
intensity that was achieved by searching for the brightest speckle. An average Intensity 
increase of only 3.2 was obtained (dotted line). Solid line, theoretical best case speckle 
average using adaptive optics, assuming that all geometrical aberrations are corrected 
perfectly. All measurements are normalized to the diffuse background. 



agating in the medium and demixed again at the target position to form a focus. This 
focus was a factor of 22 more intense than the typical intensity and a factor of 6 more 
intense than the maximum that was achieved by searching for the brightest speckle. 

The experiment was repeated for targets at different depths inside a 29±3 ^im thick 
sample. In Fig. 5.4 the target intensity that was reached with channel demixing is 
compared to the maximum intensity reached by translating the sample. Channel 
demixing increases the intensity at the target by a factor of on average 20.4, regardless 
of the depth of the target. It makes no difference whether the sample was meticu- 
lously placed on the brightest speckle first or whether the channel demixing proce- 
dure was started without a preceding scan of the sample position. For completeness, 
we also plot the calculated residual ballistic focusing component that can be expected 
when the sphere is placed exactly in the geometrical focus. For this calculation, we 
optimistically assumed that all geometric aberrations are perfectly corrected by an 
adaptive optics system[22]. As can be seen in Fig. 5.4, adaptive optics could be effec- 
tive for targets that are only a few mean free paths deep. When the target lies deeper, 
scattering on the sub -wavelength structure of the material destroys even a perfectly 
corrected geometrical focus. Using channel demixing, however, we consequently find 
that the increase does not depend on the position of the target. Therefore, we con- 
clude that scattering in no way limits the depth at which light can be delivered inside 
the metamaterial. 

Theoretically, in the current experimental configuration, channel demixing has the 
possibility to achieve an intensity ratio that is approximately four times as high as 
the measured value of 20.4 (see Section 5.B). In practice, however, we are limited by 
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photobleaching of the probe. We expect to achieve even higher intensity ratios as 
more stable probes, such as quantum dot assemblies [29], become available. 

The measured set of complex propagation amplitudes characterizes propagation 
of light towards a sub-wavelength target inside the metamaterial. For this measure- 
ment no physical access to the target is required; it is not even required to know 
where the target probe exactly is. In both ordered and disordered metamaterials these 
amplitudes contain information about the structure and the optical quality of the 
metamaterial. [20, 30, 31, 32] This information, which was not experimentally acces- 
sible before, makes channel demixing a unique method for characterizing any meta- 
material. Our experimental results call for the development of theoretical approaches 
that treat light propagation inside disordered photonic metamaterials. [33] 

We have shown that light can be focused on a target inside a photonic metama- 
terial that completely scrambles the incident waves. Since the illumination is selec- 
tive, channel demixing can be used to increase the contrast with a fluorescent back- 
ground. A different application would be to selectively activate photochemical pro- 
cesses by concentrating light at a point inside a scattering material. Finally, the use 
of a fluorescent probe permits selectively focusing light on fluorescent areas, which 
is highly desirable in biological imaging, and may be a key technology for selectively 
illuminating individual cells in disordered biomaterials such as human skin tissue. 



5.A Experimental details 

5.A.1 Apparatus 

The goal of the experiment is to illuminate a disordered medium in such a way that 
light focuses at a predefined target deep inside. The experimental setup consists of 
a spatial light modulator that is used to adjust the illumination, and a confocal fluo- 
rescence microscope that indirecfly measures the intensity at the target. A detailed 
schematic of the setup is given in Fig. 5.5. Light from a continuous wave 532 nm 
laser (Coherent Compass M315-100) is expanded using a 30x beam expander. The 
expanded beam illuminates a Holoeye LC-R 2500 spatial light modulator (SLM) with 
a resolution of 1024 x 768 pixels. The pixels are grouped into square segments of 20 
X 20 pixels. We have independent control over the phase and amplitude of each of 
the segments. [34] During this experiment, the amplitude of all segments is kept con- 
stant. The SLM is imaged onto the back aperture of a microscope objective using a 1:2 
demagnifying telescope. In total, 640 segments fall within the aperture of the objec- 
tive, which means that we control 640 degrees of freedom. The microscope objective 
(Zeiss infinity corrected 63x, NA=0.95 Achroplan air) focuses the light on the sample, 
thereby mapping the segments to channels in k-space. The sample is mounted on 
a XYZ-nanopositioning stage (Physik Instrumente P-611.3S NanoCube) to allow for a 
precise 3D-scan of the geometrical focus. 

The returned fluorescence light is separated from the excitation light using a com- 
bination of a dichroic mirror (Semrock FF562-Di02-25x36) and a bandpass filter 
(Semrock FFOl-593/40-25). An electron multiplying CCD camera (Andor Luca DL- 
658M) images the fluorescence emission. A computer sums the intensity in the cam- 
era image to find the total power in the diffuse spot of fluorescence light. For each 
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Figure 5.5: Experimental setup. Thie left part consists of a laser (532 nm) thiat Is expanded by a 30x 
beam expander and modulated with a spatial light modulator (SLM). A 1 :2 demagnifying 
telescope images the SLM on the back aperture of a 63x microscope objective. The ob- 
jective focuses the light on a sample that is mounted on a XYZ piezo positioning stage. A 
dichrolc mirror (DM) and a bandpass filter (F) block the excitation light. The fluorescence 
emission Is Imaged on an electron multiplying CCD camera (EMCCD). D, Iris diaphragm; 
PBS, polarizing beam splitter cube; Li, Lz, L3, lenses with a focal distance of 200 mm, 100 
mm and 150 mm, respectively. Some folding mirrors and beam attenuation optics were 
omitted. 



of the segments, measurements are performed at phase delays of 0, 2n/5, An/ 5, 6n/5, 
and 871/5 and a sinusoid is fitted to these measurements. Before advancing to the 
next segment, the phase delay is reset to 0. 



5.A.2 Measurement sequence 

Channel demixing is achieved by measuring the propagation coefficients from the 
incident channels towards the embedded target. After all propagation coefficients 
are measured, it is known which linear combination of scattering channels focuses at 
the target. To quantify the effectiveness of this method, we measure the intensity at 
the probe position before and after the demixing procedure. 

Since the field inside the medium is a completely random speckle pattern, the ini- 
tial intensity is a random variable. A single measurement will not provide a reliable 
reference. Therefore, we measure the speckle averaged initial intensity Iq by gener- 
ating 500 different random combinations of incident channels and determining the 
average fluorescence intensity. This way, fluctuations caused by volume speckle are 
averaged out and a reliable measure for the initial intensity is found. 

Next, the propagation constants for the different scattering channels are measured. 
This measurement takes approximately 18 minutes. During these measurements, the 
dye in the sphere photobleaches, which causes the response of the probe to decrease. 
Since we are interested in how much excitation light is focused on the target position. 
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Figure 5.6: Evolution of tine fluorescence intensity during the measurements. The measurement pro- 
cedure has three phases. 0-130 s, reference measurement; 130-1187 s, measurement of 
propagation coefficients; at 1 1 87 s, construction of optimal wavefront; 1 1 87-2244 s, fluo- 
rescence measurement with optimal wavefront. 



we have to take into account the reduced response of the probe particle 
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(5.1) 



where /opt is the fluorescence intensity directly after constructing the optimal com- 
bination of scattering channels. 7i is the fluorescence intensity at the start of this 
measurement step (directly after the reference measurement), and I2 is the intensity 
at the end of the measurements (directly before the SLM is programmed with the op- 
timal configuration). The ratio h/Iz - which is between 1.0 and 1.8 when the sample 
is positioned at a typical position and between 1.5 and 3.9 when the sample is posi- 
tioned at the brightest speckle - corrects for the fact that photobleaching decreases 
the sensitivity of the probe. 

Figure 5.6 shows how the fluorescence intensity evolves during the measurement 
sequence. During the reference measurement, the intensity fluctuates rapidly. In the 
next step, the propagation constants are measured. Since we only cycle the phase of 
a single channel at a time, the illuminating wavefront is practically constant. During 
this step the fluorescence intensity decreases slowly due to photobleaching. After the 
measurements have been performed for all segments, the phase modulator is pro- 
grammed with the optimal phases and the fluorescence intensity jumps to a higher 
value as the incident light now interferes constructively at the target. After this jump, 
the wavefront is kept constant. The reason that the fluorescence intensity decreases is 
that the higher incident intensity on the sphere causes the dye to bleach more rapidly. 



5.A.3 Sample preparation and characterization 

We demonstrated how interference of scattered light was used to focus light deep in- 
side disordered metamaterials. Our samples are strongly scattering, that is, the mean 
free path for light is of the order of the wavelength. Inside the sample, fluorescent 
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Figure 5.7: Inverse total transmission through layers of ZnO pigment. The solid line is a fit of Eq. (5.2) 
with i = 0.72 |jm, Zei = 1.15 |im, and z^i = 1.24 |im. Error bars indicate the standard 
deviation of the measured thickness and total transmission at different positions on the 
sample. 



probes mark the target positions. Typically, the probes are over 10 mean free paths 
deep, which makes it impossible to focus light on them using conventional means 
(Also see Section 5.B.4). In this section we address the fabrication and characteriza- 
tion of these samples. 

Sample preparation 

Our samples consist of a layer of Zinc Oxide pigments (Aldrich Zinc Oxide < 1 p,m, 
99.9%). The average particle diameter is 200 nm (measured using scanning electron 
microscopy). The fluorescent probes are made from polystyrene and have a speci- 
fied diameter of 300 ± 15 nm (Duke Scientific red fluorescent nanospheres dyed with 
FireFly^**). An aqueous suspension of pigment and spheres was sprayed on a stan- 
dard microscope cover slide using a commercial airbrush (Harder & Steenbeck Evolu- 
tion). The thickness of the dried samples was measured with an optical microscope 
by making a scratch in the surface. The background fluorescence of the ZnO was 
shown to be negligible compared to the fluorescence of the dyed nanospheres. 



Measurement of the mean free path 

Transport of light in scattering media is quantified by the transport mean free path ^. 
A method for obtaining i is to measure the total (diffuse) transmission of light using 
an integrating sphere. I then follows from the relation[35] 



1 

r 



L-\- Zg\ + Zg 

i+Zel 



(5.2) 



where T is the total intensity transmission, L is the sample thickness and z^i, Zez 
are the extrapolation lengths that describe boundary effects on the front and back 
surfaces of the sample respectively. The extrapolation lengths follow from the ef- 
fective refractive index of the medium, neff-[36] We measured the total transmission 
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through a set of similar samples with varying thickness. The inverse total transmis- 
sion for an ensemble of 14 samples is plotted in Fig. 5.7. The samples range in thick- 
ness from 6^xm to 32tim. As expected, the inverse total transmission increases lin- 
early with sample thickness. We obtain a good linear fit for £ — 0.7 ± 0.2 jxm and 
rigg — 1.4 ± 0.1. The extrapolation lengths can then be calculated from rigg and we 
findZei = 1.15±0.2|jLmandZe2 = 1.24±0.2|jLm. 



Determination of tlie deptli of the fluorescent probe 

For each experiment we select a different probe sphere at a different depth inside a 
sample. To determine the depth of the probe, we image the diffuse spot of emitted 
light and fit this data with diffusion theory. The accuracy of methods that rely on in- 
verting the diffusion equation is limited to about one mean free path due to noise. [37] 
The intensity inside a disordered medium with a source at depth Zo follows from 
diffusion theory. We solved the diffusion equation for a slab of thickness L using the 
Dirichlet boundary conditions /(— Zgi) — 0, and I{L + Zez) — [38, 39]. The extrap- 
olation lengths Ze\ and Ze2 account for the reflections at the front and the back of 
the sample, respectively. A closed form solution is found using Fourier transformed 
traversal coordinates qi = {q±x,<liy)- 



idiqi.z)-- 



' r sinh(<7i [Le-z- z^i]) sinh(t?i [zp + Zei]) 
^'" D^isinh(^iL,) ^>^° 

smh[qx_[Le-Zo-Zei])smh[qi_[z + Zei]) ^ ' ' 

^'" D^isinh(^iL,) ^-^^ 



where Id is the diffuse energy density at depth z, q]_ = Jq^^ + q'l , and /in is the 
power of the source in Watts. D is the diffusion constant and Le = L + Zei + Zei- 

To obtain the depth of the source, we Fourier transform the microscope image of 
the fluorescence light and fit Eq. (5.3) with z=0 using Zo and an intensity prefactor as 
fit parameters. 



5.A.4 3-dimensional scan results 

Conventional focusing does not work in strongly scattering materials such as our 
samples. The absence of a ballistic focus was confirmed by performing an extensive 
3-dimensional scan of the sample with respect to the microscope objective. The flu- 
orescence emission strongly fluctuates with the sample position. These fluctuations 
are the result of a volume speckle field in the medium that is formed by scattered 
monochromatic light. As the sample moves, the volume speckle fluctuates randomly. 
For each depth of the geometrical focus, we record the brightest speckle intensity. 
Figure 5.8 shows the results of a typical 3-D scan. The depth of the geometrical fo- 
cus was estimated to be n^sz, with n^s — 1-4 as measured, and z the position of the 
translation stage. No sign of a ballistic focus was observed. In other words, it is not 
possible to focus light on the target using conventional means. 
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Figure 5.8: Determination of the absence of a geometrical focus. The sample is moved with respect 
to the objective (inset). For each depth, the sample is scanned in a plane parallel to the 
surface to find the maximum fluorescence intensity (solid line). In homogeneous media, 
there would be a sharp peak at the point where the fluorescent sphere is in the geometrical 
focus of the microscope objective (dotted line). In disordered photonic media, a geometric 
focus cannot be formed. The highest intensity is found when 'focusing' just below the 
sample surface. 



5.B Analysis of the channel demixing method 

The optical field £ at a point rt inside the medium is given by 



Eitb]^ d^r„G(rb,rJS(rJ 



(5.4) 



where G is Green's function for propagating from the sources S{ra) to point r^. In a 
disordered medium, G is stochastic and, therefore, initially unknown. We treat each 
segment of the phase modulator as a sheet source with amplitude A and phase 4>- 
Since the phase modulator is in good approximation illuminated homogeneously, all 
amplitudes A are equal. By integrating over the surface area S of each of the N seg- 
ments, Eq. (5.4) is discretized. 



£(r.) = Xff 



d^raGirb,ra)Ae'^", 



i^ttae'*^. 



(5.5) 
(5.6) 



When the phase of a single segment a of the phase modulator is changed, the inten- 
sity at point Tb responds as 



7(rb) = |£(rfc)p = 7ob +2^Re(£* ffc„e'^«), 



(5.7) 
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with 

lob^lEbaf+A'ltbaf, (5.8) 

JV 

Em^A^ tba'e''''"' « Eirb). (5.9) 

Since the number of segments N is large, Ebd "^ E{rb) and, therefore, equal for each 
of the segments. By repeating this procedure for all segments a, we measure the co- 
efficients tba up to an unknown common prefactor E{rb)- 

Once all coefficients tba are known, our channel dembdng method maximizes E{rb) 
by setting 0^ — —aig{tba) for all segments. This configuration of phases gives the 
global maximum of Eq. (5.6), 

N 

EmUn)^A^\tba\. (5.10) 

a 

A comparison with Eq. (5.6) shows that now all amplitudes are summed with the same 
phase and that all incident channels will interfere constructively at the target. 

5.B.1 Maximum enhancement - scalar waves, simplified case 

The theoretical increase in intensity that can be achieved with channel demixing de- 
pends on the number of modulator segments N and the statistical properties of coef- 
ficients tba ■ In disordered systems it is often allowed to assume that the coefficients 
are independently drawn from a single circular Gaussian distribution. [12, 13, 40, 41, 
42] Under this assumption, the expected ratio between the original diffuse intensity 
and the maximally achievable intensity equals 

rj = -(Af-l) + l. (5.11) 

where N can never exceed the number of mesoscopic scattering channels. The num- 
ber of contributing scattering channels depends on the depth of the target. The 
deeper the target is embedded in the scattering material, the more channels will 
contribute to the mixing process. We estimate the number of channels that con- 
tribute to the target field by examining the size of the diffuse spot of the fluores- 
cence light. The number of scattering channels in a spot with surface area A is given 
by[43, 44] Afmax = 2nA/X^. In Fig. 5.3b, the full width at half maximum of the diffuse 
spot is 15 [xm, which leads to A = 177 ^xm^ and A/max = 3.9 • 10^. In our experiment 
N = 640 ^ 3.9 • 10^ and neighboring segments were uncorrelated. 

5.B.2 Maximum enhancement - finite size probe 

Ideally, the sphere that is used as a probe is only sensitive to the field of a single po- 
larization at a single point. In reality, however, the fluorescent sphere probes the 
summed intensity of three polarizations in a finite volume. Therefore, it is as if we 
have multiple targets for which we try to increase the intensity simultaneously. As a 
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result, the extra intensity is divided over all targets. The ratio between the original 
diffuse intensity and the maximally achievable intensity now equals 

n N-l 

where M is the number of targets. We determine M by calculating how many orthog- 
onal fields can exist in the volume of the probe sphere. Using Mie theory we find 
that for spheres of 300 ± 15 nm the three electric dipole modes and the three mag- 
netic dipole modes all give an equally strong contribution to the fluorescence power. 
The sphere is too small to support higher order (quadrupole) modes, so we estimate 
that the probe contains 6 optical modes. Using Eq. (5.12) with M = 6 and N — 640, 
we find an expected theoretical enhancement of 85. The experimental enhancement 
was approximately 25% of this value. We believe that the most important reason for 
not reaching the theoretical maximum enhancement is photobleaching. The chan- 
nel demixing method favors the mode that initially is the brightest. Unfortunately, 
this target mode also photobleaches the fastest. Since we correct only for the average 
photo bleaching rate (see Section 5.A.2), the measured enhancement will systemati- 
cally be lower than the actual enhancement. 

5.B.3 Maximum enhancement of speckle scan 

When the sample is moved, the speckle pattern in the medium fluctuates. As a result, 
the probe fluorescence will fluctuate too. Therefore, when a 3D-scan of the sample 
position is performed, at some points the intensity will be somewhat higher than the 
diffuse background. We used numerical simulations to calculate the expected maxi- 
mum enhancement that this scanning technique can achieve. We assumed Rayleigh 
statistics for the speckle field and calculated the summed intensity in 6 independent 
optical modes at Nr different positions. The simulation was repeated 1000 times and 
for each run the maximum speckle intensity was stored. The number of independent 
speckles, N^, is limited by the diffuse intensity envelope of the light inside the mate- 
rial. We estimate that Nr is between 1000 and 10000, depending on the depth of the 
probe. With these numbers, the simulations predict that the brightest speckle is 2.9 
to 3.4 times as intense as the diffuse background. 

5.B.4 Diffuse and ballistic intensities inside the medium 

Diffuse intensity 

When a slab of diffusive material is illuminated, the incident light can be described 
by a diffuse source at depth Zo=lm the medium. [45] The maximum diffuse intensity 
at depth z is found by inverse Fourier transforming Eq. (5.3) with x = y = 0. 

1 r 

/rf(x = 0,y = 0,z)=— dqj_Idiqi,z)Siq±)q±. (5.13) 

2^ Jo 

Where S{q±) accounts for the spatial extent of the source. We modelled the source 
as a Gaussian spot that encompasses 640 scattering channels; S{q±] = exp(— ^|/0.43). 
Eq. (5.13) was evaluated numerically. The diffuse intensity was used to normalize the 
theoretical curve for the perfect ballistic intensity in Fig. 5.4. 
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Perfect ballistic intensity 

We assume that all geometrical abberations are perfectly corrected and that all the 
incident light is directed towards a sphere in the geometrical focus. While travelling 
through the disordered medium, scattering and diffraction on sub -wavelength struc- 
tures attenuate the incident beam according to the Lambert-Beer law, 

4 = — , (5.14) 

Av 

where h is the ballistic energy density at the sphere, A is the cross-section area of 
the sphere, z is the depth of the sphere, v is the velocity of light inside the medium 
and tsc is the scattering mean free path for light. The parameters tsc and v have not 
been measured directly for our samples. We can, however, give an upper limit for 
the ballistic intensity by overestimating Isc — i and underestimating v — 3D/1. The 
ballistic intensity was normalized by the diffuse intensity and plotted in Fig. 5.4. 
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Chapter 6 

Exploiting the potential of 
disorder in optical communication 



Scattering of light is a long recognized problem that causes losses in optical telecom- 
munication, and decreases the resolution of imaging systems[l]. However, impres- 
sive experiments utilizing ultrasound and radio waves have shown that, for these 
waves, scattering and interference can be used as an advantage. [2] Disordered scat- 
tering is used for ultrafast radiographic communication [3] and for sub-wavelength 
focusing of microwaves [4]. 

Unfortunately, these revolutionary methods have so far not been available in the 
optical regime. The reason is that recording and playing back electrical fields at op- 
tical frequencies is close to impossible. Therefore, time-reversal - which is highly ef- 
fective for sound and radio waves - is not a realistic option. In this chapter, we show 
that the potential of scattering can be exploited without time-reversal. A self-learning 
optical system uses scattering to improve its resolution beyond the limit imposed by 
diffraction at its aperture. Thanks to the increased resolution, multiple communica- 
tion channels can coexist in the space that is normally needed to support one unscat- 
tered communication channel. 

We first explain the principle of the experiment and present our results. Then, in 
Section 6.2, we give a more detailed description of the apparatus and the experimen- 
tal procedure. 



6.1 Increasing the information density 

Figure 6.1a shows an optical system where a lens focuses monochromatic light on a 
camera. Due to the limited numerical aperture (NA=l/200) of the system, the spot 
will at best have a diffraction limited width of 60jim (full width at half maximum). 
While aberrations in the optics can efficiently be removed using adaptive optics [5], 
the resolution of an adaptive optics system is still limited by its numerical aperture. 

We choose three individual pixels of the camera to define three independent re- 
ceivers. The spacing between the pixels is less than the diffraction limit of the lens. 
When the beam path is blocked by a disordered scattering medium (a 9 ^xm- thick 
layer of airbrushed white paint), the image on the camera changes dramatically. In- 
stead of a focal spot, we now record a disordered speckle pattern that is the result 
of random interference. Although the scattered light does not focus, its features are 
much sharper than the size of the original diffraction limited spot because the scat- 
tering process increases the effective numerical aperture of the optical system. 

To control the propagation of light through the disordered environment, we devel- 
oped a learning algorithm that acquires an optical fingerprint (a spatial amplitude 
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Figure 6.1: a) An optical system with a limited numerical aperture cannot focus light sharply, all three 
detectors receive the same signal, b) Scattering improves the diffraction limit. The incident 
light is spatially shaped to match a scattering fingerprint. Scattered light focuses sharply 
to a point. Now, each detector can be addressed individually c) Measured transmitted 
intensity for eight different linear combinations of incident wavefronts. Any three-bit pattern 
can be generated in the area previously needed for one bit. 
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and phase pattern) for each of the three receivers. The optical fingerprint uniquely 
describes propagation through the scattering medium towards the receiver. When 
the wavefront of the incident light matches this fingerprint, the scattered light con- 
verges to a focus at the receiver (Fig. 6.1b) and simultaneously the intensity on the 
other two receivers decreases. 

We used a recently developed technique to spatially shape both phase and am- 
plitude of a laser beam with a high accuracy. [6] This technique allowed us to con- 
struct any linear superpositions of fingerprints. Propagation through the scattering 
medium separates the superposition into separated foci on the receivers. All patterns 
of three bits can be generated, as is shown in Fig. 6.1c. A thin layer of white paint has 
increased the information density by a factor of three. 

In conclusion, our method for matching optical fingerprints takes highly successful 
ideas from microwave technology and introduces them to the optical regime. Apart 
from improving optical communication, applications include sub -wavelength focus- 
ing with ordered [7] or disordered[4] metamaterials. 

6.2 Experimental details 

The setup consists of a low numerical aperture optical system and a wavefront syn- 
thesizer that spatially modulates light from a 632.8 nm Helium-Neon laser (see Fig- 
ure 6.2). Wavefront modulation is achieved using a Holoeye LC-R 2500 spatial light 
modulator (SLM). The SLM is of the twisted nematic type and by itself does not allow 
independent control over the phase and amplitude of the reflected light. We devel- 
oped a technique that combines four pixels into a macropixel of which both phase 
and amplitude can be controlled fully. [6] 

Without scattering, the optical system is not able to create a sharp focus, regardless 
of the wavefront that the wavefront synthesizer produces. When we place a layer of 
strongly scattering paint in the optical path, light is scattered. The scattered light 
forms a random speckle pattern with sharp features. To focus light to a single point, 
we generate a wavefront that precisely matches the scattering matrix of the layer of 
paint. The scattered light now interferes constructively in a single point. 

To construct the matching wavefront, or 'fingerprint', we use a method that is simi- 
lar to the technique described in Chapter 3. One by one the phase of a segment of the 
phase modulator is cycled from to 2n. Meanwhile, the intensity at the desired focal 
point is monitored with a camera. From the resulting interference signal we deduce 
the phase and amplitude of that segment's contribution to the target focus. This mea- 
surement is repeated for all segments. From these measurements, the unique [8] opti- 
cal fingerprint that causes constructive interference in the target focus is constructed. 
In contrast to the method presented in Chapter 3, we also control the amplitude of the 
incident light. We set the amplitude of each of the segments to be proportional to the 
measured amplitude of that segment's contribution to the target focus. In this way the 
overlap of the incident field with the fingerprint of the medium is optimal. Segments 
that contribute little to the focus are effectively switched off so that the background 
intensity outside of the focus is reduced and contrast is improved. Theoretically, am- 
plitude and phase modulation increases the contrast with the diffuse background by 
a factor of 4/n with respect to phase only modulation. Moreover, thanks to the com- 
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Figure 6.2: Schematic of thie experimental setup. The first part of the setup consists of a wavefront 
synthesizer that modulates amplitude and phase of a 632.8 nm Helium-Neon laser. The 
laser beam is expanded by a 30x beam expander and modulated with a spatial light mod- 
ulator (SLM). A 1 :5 demagnifying telescope images the SLM on the surface of lens L3. A 
pinhole in the focal point of the telescope (Ai) acts as a spatial filter, which allows the SLM 
to modulate both amplitude and phase (see text). The second part is the actual experiment, 
consisting of a camera (CCD) that is placed in the focal plane of lens L3. The numerical 
aperture of the lens is lowered using a 2.1 -mm diameter aperture (A2). At 26 mm from the 
camera a 9|im-thick strongly scattering layer of white paint is placed. Aj, A2, apertures; 
PBS, polarizing beam splitter cube; Li, L2, L3, lenses with a focal distance of 250 mm, 
50 mm and 200 mm, respectively. Some folding mirrors and beam attenuation optics were 
omitted. 



bined amplitude and phase control, fingerprints for different foci can be superposed 
to generate any desired bit pattern. 
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Chapter 7 

Phase control algorithms for 
focusing light through turbid 
media 



Abstract: Light propagation in materials witli microscopic inliomogeneities is af- 
fected by scattering. In scattering materials, such as powders, disordered metamate- 
rials or biological tissue, multiple scattering on sub-wavelength particles makes light 
diffuse. Recently, we showed that it is possible to construct a wavefront that focuses 
through a solid, strongly scattering object. The focusing wavefront uniquely matches 
a certain configuration of the particles in the medium. To focus light through a tur- 
bid liquid or living tissue, it is necessary to dynamically adjust the wavefront as the 
particles in the medium move. Here we present three algorithms for constructing a 
wavefront that focuses through a scattering medium. We analyze the dynamic be- 
havior of these algorithms and compare their sensitivity to measurement noise. The 
algorithms are compared both experimentally and using numerical simulations. The 
results are in good agreement with an intuitive model, which may be used to develop 
dynamic diffusion compensators with applications in, for example, light delivery in 
human tissue. 
[This chapter has been accepted for publication in Optics Communications.] 



Materials such as paper, white paint or human tissue are non-transparent because 
of multiple scattering of light. [1, 2, 3] Light propagating in such materials is diffuse. 
Recently, we have shown that coherent light can be focused through diffusive me- 
dia yielding a sharp, intense focus. [4] Starting with the situation where a scattering 
object (a layer of Ti02 pigment with a thickness of approximately 20 transport mean 
free paths) completely destroys the spatial coherence of the incident light (Fig. 7.1a, 
7.1c), we controlled the incident wavefront to exactly match scattering in the sample. 
Afterwards, the transmitted light converged to a tight, high contrast focus (Fig. 7.1b, 
7. Id). These matched wavefronts experience inverse diffusion, that is, they gain spa- 
tial coherence by travelling through a disordered medium. 

For a given sample of scattering material, there is a unique incident wavefront that 
makes the object optimally focus light to a given point. Like a speckle pattern, this 
wavefront is disordered on the scale of the wavelength of light. This wavefront cannot 
be constructed from a small number of smooth base functions, which unfortunately 
renders the efficient algorithms used in adaptive optics (see e.g. [5]) ineffective. In 
Ref. (author?) 4, we presented an algorithm that finds the optimal wavefront when 
the sample is perfectly stationary and the noise level is negligible. To find applications 
in, for example, fluorescence excitation or photodynamic therapy, the wavefront has 
to be adjusted dynamically as the scatterers in the sample move. In this paper, we 
present two additional algorithms, that dynamically adjust the wavefront to follow 
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Figure 7.1 : Principle and experimental results of inverse wave diffusion, a) A multiply scattering object 
destroys the spatial cofierence of incident light, b) When the same object is illuminated with 
a specially constructed matching wavefront, the transmitted light focuses to a tight spot, c) 
Recorded intensity transmission of an unshaped wave through a 10|im thick layer of Ti02 
pigment, d) Intensity transmission through the same sample with a shaped wavefront. 



changes in the sample. The performance of the algorithms is in good agreement with 
numerical simulations and with an analytical model. We show that the new algo- 
rithms are superior to the original algorithm when the scatterers in the sample move 
or when the initial signal to noise ratio is poor. 

Wave diffusion is a v«dely encountered physical phenomenon. The use of mul- 
tiple scattered waves is the subject of intensive study in the fields of, for instance, 
ultrasound imaging[6, 7, 8], radio and microwave antennas[9, 10], seismography[ll], 
submarine communication [12], and surface plasmons[13]. While the algorithms dis- 
cussed in this paper were developed for spatial phase shaping of light, they can be 
used for any type of wave and apply to spatial phase shaping as well as to frequency 
domain phase shaping (also known as coherent control, see e.g. [14]) as the concepts 
are the same. 

This article is organized as follows. First the key concepts of inverse diffusion are 
introduced and the three different algorithms are presented. Then the experimen- 
tal apparatus is explained and the measured typical performance of the algorithms is 
compared. In the subsequent section, we will compare the experimental results with 
numerical simulations and analyze the data in terms of noise and stability of the scat- 
terers. Finally, we will analytically explain the characteristic features of the different 
algorithms and discuss their sensitivity to noise. 
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Figure 7.2: Feedback loop for achieving inverse diffusion. An incident monochromatic beam is shaped 
using a spatial light modulator (1) and projected on a non-transparent multiply scattering 
object (2). A detector (3) detects the amount of transmitted light that reaches the target 
area. A feedback algorithm (4) uses the signal from the detector to program the phase 
modulator. Before the algorithm is started, the transmitted light forms a random speckle 
pattern. The algorithm changes the incident wave to increase the intensity in the target 
area. After a few iterations, the transmitted light focuses on the target. 



7.1 Algorithms for inverse diffusion 

The key elements of an inverse diffusion setup are a multiply scattering sample, a spa- 
tial light modulator, an optimization algorithm and a detector, as shown in Fig. 7.2. 
The sample can be anything that scatters light without absorbing it. We will con- 
sider only samples that are thicker than approximately 6 transport mean free paths 
for light. Light transmitted through these samples is completely diffuse and the trans- 
mitted wavefront is completely scrambled, i.e., it has no correlation with the incident 
wavefront. [15] We assumed that the intensity is low so that non-linear effects are neg- 
ligible. 

The incident wavefront is constructed using a spatial phase modulator. The mod- 
ulator consists of a 2D -array of pixels that are grouped into N equally sized square 
segments. A computer sets the phase retardation for each of the segments individ- 
ually to a value between and 2?!. The optimization algorithms program the phase 
modulator based on the detector output. Initially, each of the segments randomly 
contributes to the field at the detector. The algorithms find the unique configura- 
tion of retardations for which all contributions are in phase at the target. As follows 
from the triangle inequality, the magnitude of the field at the detector is at a global 
maximum for this optimal configuration. Since the sample completely scrambles the 
incident wavefront, all segments of the wavefront are scattered independently and 
the optimal wavefront will not be smooth. 

Behind the sample is a detector that provides feedback for the algorithm. The de- 
tector defines the target area where the intensity is maximized. The field at the detec- 
tor is the result of interference from scattered light originating from the different seg- 
ments of incident wavefront. When the phase of one or more segments is changed, 
the target intensity responds sinusoidally. We sample the sine wave by taking 10 mea- 
surements. The process of capturing a single sine wave and possibly adjusting the 
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Figure 7.3: Principle used in tine three different optimization algorithms, a) For the stepwise sequential 
algorithm, all segments are addressed sequentially (marked squares). After the optimal 
phase is measured for all segments, the modulator is updated to construct the optimal 
wavefront (light gray squares), b) The continuous sequential algorithm is equal to the first 
algorithm, except that the modulator is updated after each iteration, c) The partitioning 
algorithm randomly selects half of the segments and adjusts their overall phase. The mod- 
ulator is updated after each measurement. 



phase modulator accordingly is called an iteration. 

The amount of control we have over the propagation of light in the disordered sys- 
tem is quantified by the signal enhancement. The enhancement 17 is defined as 



17 = 



Uo)' 



(7.1) 



where I^ is the intensity in the target after optimization and (Iq) is the ensemble av- 
eraged transmitted intensity before optimization. In a perfectly stable system, the 
enhancement is proportional to N [4], meaning that the more individual segments 
are used to shape the incident wavefront, the more light is directed to the target. In 
practice, however, the enhancement is limited by the number of iterations that can 
be performed before the sample changes too much. We define the persistence time 
Tp as the decay time of the field autocorrelate of the transmitted speckle, which is 
a measure for the temporEil stability of the sample. The persistence time depends 
on the type of sample and on environmental conditions. Typical values of Tp range 
from a few milliseconds in living tissue[16] to hours for solid samples in laboratory 
conditions. The other relevant timescale is the time required for performing a single 
iteration, Tt . In our experiments, we operate the phase modulator at just below 10 Hz 
and take ten measurements for each iteration; we have i; ^ 1.2s. 

We will now present three algorithms we used to invert wave diffusion. The advan- 
tages and disadvantages of the algorithms are discussed briefly and will be analyzed 
in detail later. 



7.1 .1 The stepwise sequential algorithm 

The stepwise sequential algorithm that was used in Ref. [4] is very straightforward. It 
is based on the fact that the field at the detector is a linear superposition of the contri- 
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butions from all segments. This means that we can construct the optimal wavefront 
by optimizing each of the segments individually. The computer consecutively cycles 
the phase of each of the N segments from to 2n. The feedback signal is monitored 
and the phase for which the target intensity is maximal is stored. At this optimal 
phase, the contribution of this single segment is in phase with the background field 
(the mean contribution coming from all other segments). The phase retardation of 
the segment is reset to before continuing with the next segment. This way, the back- 
ground field remains unchanged. Only after all iterations are performed, the phase of 
each segment is set to this optimal value (see Fig. 7.3a). Now, all contributions have 
the same phase as the original background field. 

In absence of measurement noise or temporal instability, algorithm 1 is guaran- 
teed to find the global maximum in the lowest number of iterations possible. How- 
ever, when NTi "> Tp, the speckle pattern decorrelates before all measurements are 
performed and the algorithm will not work. Therefore, it is important to adjust the 
number of segments to the persistence time. 



7.1.2 The continuous sequential algorithm 

The continuous sequential algorithm is very similar to the stepwise sequential algo- 
rithm except for the fact that the phase of each segment is set to its maximum value 
directly after each measurement (see Fig. 7.3b). This approach has two advantages. 
First of all, the algorithm runs continuously and dynamically follows changes in the 
sample's scattering behavior. Furthermore, the target signal starts to increase directly, 
which increases the signal to noise ratio of successive measurements. It is still neces- 
sary to adjust N to the persistence time Tp. 



7.1.3 The partitioning algorithm 

The two sequential algorithms change the phase of one segment at a time. As one 
segment contributes only a small fraction of the wavefront, the signal to noise ratio is 
low and the initial convergence is slow. As an alternative, we propose a partitioning 
algorithm that simultaneously changes the phase of a whole group of segments (see 
Fig. 7.3c). The amplitude of the interference signal is maximal when the group con- 
tains half of the segments. At each iteration the segments are divided randomly over 
two equally sized partitions. Then, the target intensity is maximized by changing the 
phase of one partition with respect to the other. This process is repeated indefinitely, 
and the modulator is randomly repartitioned every time. 

Since the phase of many segments is changed simultaneously, the initial increase 
in intensity will be fast and the feedback signal will be high. Therefore, this algorithm 
is expected to be less sensitive to noise and to recover from disturbances more rapidly. 



7.2 Experiment 

The different algorithms were tested experimentally using the setup shown in Fig. 7.4. 
In our case the scattering medium is a 10 ^xm thick layer of rutile TiOz pigment [17] 
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Figure 7.4: Experimental apparatus used for inverting diffusion. Light from a HeNe laser is spatially 
modulated by a liquid crystal spatial light modulator (SLM). Wave plates and a polarizer 
(P) are used to generate and select the polarization state for which the modulator works 
in phase-mostly mode. The shaped beam is demagnified by a factor of 3 (telescope not 
shown) and focused on the sample (S). A reference detector monitors the total intensity 
falling on the sample. A microscope objective, a polarizer and a CCD-camera are used to 
detect the intensity in the target focus, a few millimeters behind the sample. M, mirror; BS, 
non-polarizing 50% beam splitter. 



with a mean free path of 0.55 ± 0.1 ^Lin, determined by measuring the total transmis- 
sion at a wavelength of 632.8 nm. This sample is illuminated by a 632.8 nm HeNe 
laser. The laser beam is expanded and spatially modulated by a Holoeye LC-R 2500 
liquid crystal light modulator (LCD). The modulator is a reflective twisted nematic 
liquid crystal on silicon (TN-LCoS) device. Such devices are designed for intensity 
modulation. However, by choosing an appropriate combination of the incident po- 
larization and the analyzed polarization, the LCD can be used as a phase modulator 
with a minimal residual intensity modulation[18]. We characterized the LCD using 
the method described in [19] and select an optimal combination of elliptical polar- 
izations to achieve a phase-mostly modulation mode with 2n phase modulation and 
a maximum of 21% residual intensity modulation. The shaped beam is focused on 
the sample using a 63x objective with a numerical aperture (NA) of 0.85. A 20x objec- 
tive (NA=0.5) images a point that is approximately 3.5 mm behind the sample onto 
a 12-bit CCD camera (Allied Vision Technologies Dolphin F-145B). This point is the 
target area where we want the light to focus. A computer integrates the intensity in 
a circular area with a radius of 20 pixels (corresponding to 129 [xm in the focal plane 
of the objective). This target area is smaller than a typical speckle spot. Using this 
signal as feedback, the computer programs the phase modulator using one of the al- 
gorithms described above. 

We first run the three different algorithms with TV — 52. Since Tp/I} » 52, we do 
not expect to see decoherence effects. In total, 208 iterations were performed, which 
means that the sequential algorithms ran four times consecutively. The results of the 
optimization procedures is shown in Fig. 7.5a. Although the three algorithms reach 
the same final enhancement of intensity, there are significant differences between 
the algorithms. The enhancement for the stepwise sequential algorithm increases 
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Figure 7.5: a) Typical runs of the stepwise sequential algorithm (dotted curve), the continuous sequen- 
tial algorithm (dashed curve) and the partitioning algorithm (solid curve). All algorithms 
were run with N = 52. The sequential algorithms were repeated four times, b) Simulation 
results for N = 52 averaged over 64 runs. The simulation captures the main features of the 
three algorithms, but predicts a higher maximum enhancement. 
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Figure 7.6: a) Typical runs of the stepwise sequential algorithm (dotted curve), the continuous sequen- 
tial algorithm (dashed curve) and the partitioning algorithm (solid curve). All algorithms 
were run with N = 1804. The enhancement are plotted on a logarithmic scale, b) Simula- 
tion results for N = 1804 averaged over 84 runs. 



in discrete steps because the phase modulator is only reprogrammed every N itera- 
tions. During the first N iterations, the target signal is low and the algorithm suffers 
from noise. The saturation enhancement is reached after the second update (after 
2N iterations). The continuous sequential algorithm and the partitioning algorithm 
both start updating the wavefront immediately and, therefore, have a higher initial 
increase of the signal. The continuous sequential algorithm is the first algorithm to 
reach the saturation enhancement (after N iterations). The partitioning algorithm 
has the fastest initial increase in the target signal. It is, however, the last algorithm to 
reach the saturation enhancement since the final convergence is very slow. 

When the number of segments in the wavefront is increased, we expect to find a 
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higher target intensity. Figure 7.6a shows the experimental results for N — 1804 on a 
logarithmic scale. The final intensity enhancement is approximately 40 times higher 
than in Fig. 7.5a. A further difference is that in this situation the effects of decoher- 
ence are no longer negligible. This effect is most clearly visible with the stepwise 
sequential algorithm. The phase modulator is updated after each N iterations and 
between the updates the intensity decays exponentially with a 1/e decay of about 
Tp/Tl^ 5000 iterations. 

The convergence behavior of the three algorithms is similar to the experiment that 
is shown in Fig. 7.5. The partitioning algorithm clearly causes a higher signal en- 
hancement during the first 1000 iterations. The initial increase in the enhancement is 
linear with a slope of 0.37. Initially, this linear increase is far superior to the quadratic 
increase obtained with the continuous sequential algorithm. 

We conclude that both new algorithms are valuable improvements over the original 
stepwise sequential algorithm. These algorithms are far less sensitive to noise and 
the target signal is kept at a constant value even in the presence of decoherence. The 
partitioning algorithm has the fastest initial increase and therefore will recover from 
disturbances most rapidly. 



7.3 Simulations 

In order to obtain a better understanding of the effects of noise and fluctuations on 
the performance of the different algorithms we perform numerical simulations. The 
disordered medium is represented by a transmission matrix with elements drawn 
from a circular Gaussian distribution (more details on the matrbc representation can 
be found below). Decoherence is modelled by adding a small perturbation to the 
transmission matrix after every measurement. Finally, measurement noise is includ- 
ed by adding a random value to the simulated detector signal. 

Figures 7.5b and 7.6b show the simulated enhancement for a system with Tp/Ti — 
5000. Every iteration, ten measurements are performed for phase delays between 
and 2n. To each of these measurements, Gaussian noise with a standard deviation of 
0.3/o was added. The magnitude of the noise is comparable to experimental observa- 
tions. The three different algorithms were run with N — 52 and N — 1804 to simulate 
the experiments shown in Fig. 7.5a and Fig. 7.6a. 

The simulations are in good qualitative agreement with the experimental data. The 
result for the stepwise sequential algorithm shows that the effects of noise and deco- 
herence are simulated realistically. Furthermore, the initial signal increase and the 
long time convergence behavior correspond to the measured results. The only signif- 
icant difference is the 20% to 50% higher enhancement reached in the simulations. A 
possible explanation for this difference is the residual amplitude cross-modulation in 
our phase modulator. Due to this cross-modulation, the optimal wavefront cannot be 
generated exactly. Furthermore, the amplitude modulation decreases the accuracy of 
the measurement of the optimal phase. The partitioning algorithm is less sensitive to 
this last effect since the cross-modulation is averaged over many segments with dif- 
ferent phases. Since the simulations capture the overall behavior of the algorithms 
very well, we can use them to extrapolate to situations with a lot of noise and strong 
decoherence or, on the other hand, to perfectly stable systems. 
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7.4 Analytical expressions for the enhancement 

In this section, we analyze the performance of the algorithms with analytical theory 
and compare these results to the simulations. We describe scattering in the sample 
with the transmission matrix elements, fm„. This matrix couples the fields of the 
incident light and the transmitted light. 
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where the 0n is the phase of the nXh segment of the phase modulator. Assuming that 
the modulator is illuminated homogeneously, all incoming channels carry the same 
intensity. We write A„ — l/\/]V to normalize the total incident intensity. Elements 
Ei,E2,... correspond to single scattering channels of the transmitted light. Since we 
are interested in focusing light to a single spot, we need to consider only a single 
transmission channel, Em . The intensity transmitted into channel m is given by 
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Regardless of the values of the elements tmn of the transmission matrix, the intensity 
\Em\^ has its global maximum when the phase modulator exactly compensates the 
phase retardation in the sample for each segment, i.e. 0„ = — arg(fm„). The target 
intensities before optimization (/q) and after an ideal optimization Umax) are given by 
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For a disordered medium the elements of tmn are independent and have a Gaus- 
sian distribution. [20, 21, 22, 23] Rewriting Eq. (7.5) gives 
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where the angled brackets denote ensemble averaging over disorder. Eq. (7.7) pre- 
dicts that the expected maximum enhancement for an ideally stable, noise free sys- 
tem linearly depends on the number of segments N. For Af » 1, we have r] ^ tzN /A. 



7.4.1 Performance in fluctuating environments 

In reality, the sample will not be completely stable. Whether this instability is due 
to a drift of the sample position, movement of the scatterers, changing humidity or 
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any other cause, the transmission matrix will fluctuate over time. In the simulations, 
we modelled decoherence by repeatedly adding a small perturbation to each of the 
matrix elements. 

tmn—' I (fmn + '^)i (7.8) 

Vl+52 

where E, is drawn from a complex Gaussian distribution with mean and standard 
deviation 5. The prefactor normalizes the transformation so that (|f P) remains con- 
stant. By substituting the continuous limit of Eq. (7.8) in Eq. (7.6), we find an analytic 
expression for the effect of decoherence. 
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where r„ is the time that has past since the phase of segment n was measured. This 
simple model explains the exponential decay of the intensity that was observed in the 
measurements (see Fig. 7.6a) and predicts a decay time of Tp — Til 5'^. 

We now calculate the maximum enhancement that can be reached with the con- 
tinuous sequential algorithm in the presence of decoherence. Because the phases of 
the segments are measured sequentially, at any given time the values for r„/I; are 
equally spaced between 1 and N. From Eq. (7.9) we find a maximum intensity en- 
hancement of 



The mcLximum enhancement for both sequential algorithms is the same. However, 
since the stepwise algorithm only updates the projected wavefront after N iterations, 
the enhancement decreases exponentially between updates. 

In Fig. 7.7 the enhancement for different values of N is plotted versus Tp/Tj. When 
the persistence time is large {TplTi:> N), decoherence effects do not play a role and 
the enhancement linearly depends on N as was seen in Eq. (7.7). For Tp/Tj < N, how- 
ever, the enhancement decreases because the speckle pattern decorrelates before all 
iterations are performed and the enhancement drops to zero. As a consequence, the 
sequential algorithms only perform optimal when N is adjusted to Tp. When Tp is 
known a-priori, this optimum for N can be found by maximizing Eq. (7.10). We find 

N^pt^WTp/Ti, (7.11) 

where W ^ 2.51 is the solution of exp{W/2) — 1 + W. The maximal enhancement 
achievable with sequential algorithms follows by substituting Eq. (7.1 1) into Eq. (7.10) 
and equals rjopt — 0.640 Tp/T;. 

With the partitioning algorithm rj increases by 1/2 each iteration of the algorithm 
(see Section 7.A) . As long as N:> Tp/Tj, the enhancement saturates at rj—Tp /(2 Tj)+l, 
when the increase is exactly cancelled by the effect of decoherence. The most impor- 
tant difference with the sequential algorithms, is that the enhancement reached with 
the partitioning algorithm does not depend on N. In Fig. 7.7 it can be seen that the 
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Figure 7.7: Theoretical maximum enhancement as a function of coherence time for different algorithms. 
The solid lines represent the maximum enhancement that can be obtained using the se- 
quential algorithms. The enhancement depends on the number of segments used in the 
algorithm. The dashed line shows the enhancement for the partitioning algorithm where 

N:>TplTi. 




Figure 7.8: Simulated effect of decoherence on the sequential algorithms (solid curve) and on the 
partitioning algorithm (dashed curve). The simulations are averaged over 25 runs. Only 
when ArsaZ.SlTp/T;, the sequential algorithms perform slightly better than the partitioning 
algorithm. 



partitioning algorithm outperforms the sequential algorithms for almost all combi- 
nations of Tp and N. The sequential algorithm only give a slightly higher enhance- 
ment when they are fine-tuned for a known persistence time (JV = 2.51 Tp/T;). In most 
situations, Tp is not known a-priory or varies over time and the partitioning algorithm 
will be preferable. 

Our analytical results for all three algorithms are supported by numerical simula- 
tions (see Fig. 7.8). The simulations exactly reproduce the theoretical curves shown 
in Fig. 7.7. For the simulations of the partitioning algorithm we used N — 4096, a 
number that can easily be reached with a LCD phase modulator. Again, the partition- 
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Table 7.1 : Signal and noise characteristics of the three algorithms. The rms phase correction is a 
measure for the required sensitivity. 



ing algorithm can be seen to have good overall performance, whereas the sequential 
algorithms only work well for certain combinations of N and Tp . 

In conclusion, the maximum enhancement that can be reached linearly depends 
on the sample's persistence time. For the sequential algorithms 17 — 0.64Tp/Ti, but 
only when N is precisely adjusted to Tp . The partitioning algorithm has ri — 0.5Tp/Ti, 
as long as N is large enough. Using these analytical relations, the performance of each 
of the three algorithms in different experimental situations can easily be estimated. 



7.5 Effect of Noise 

Measurement noise affects the measured phases. Noise-induced errors in the phases 
lead to a reduction of the enhancement, rj. We will now compare the signal-to-noise 
ratio (SNR) of the three different algorithms. In a single iteration of an algorithm, the 
phase of one or more segments is varied, while the phase of the other segments is 
kept constant. The intensity at the detector equals 



/($) = /^ + /b + 2VW^cos($-$o), 



(7.12) 



where Ib is the intensity at the target originating from the modulated segments, Ia is 
the target intensity caused by light coming from the other segments, $ is the phase 
that is varied and #0 is the unknown optimal value for the phase. The last term in 
Eq. (7.12) is the signal that is relevant for measuring $0. There first two terms consti- 
tute a constant bias. 

Table 7.1 lists the magnitudes of the signal and the bias for each of the three algo- 
rithms. If the detection system is photon shot noise limited, the noise is proportional 
to the square root of the bias. When, on the other hand, constant noise sources such 
as readout noise or thermal noise are dominant, the SNR is directly proportional to 
the signal magnitude. Table 7.1 also shows the root mean square (rms) phase correc- 
tion that is applied during each iteration of the corresponding algorithm. The rms 
phase correction is a measure for the required accuracy of the measurements. 

With the stepwise sequential algorithm, Ia ^ h and on average Ib — Iq/N. Since 
the initial diffuse transmission /q is low and N can be very high, the SNR is low. The 
continuous sequential algorithm has a higher SNR since the overall intensity at the 
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detector increases while the algorithm progresses and Ia i^ tjIq. Assuming the dom- 
inant noise source is constant, the SNR will increase as the enhancement becomes 
higher. Therefore, the algorithm can be accelerated by decreasing the integration 
time of the camera as the algorithm advances. In case the detection system is photon 
shot noise limited, the SNR remains constant during the optimization since both the 
signal and the shot noise scale as ^/rjl^. 

The highest SNR is achieved with the partitioning algorithm. Since we always 
change the phase of half of the segments, Ia ^ Is ^ i7^o/2, resulting in a maximal 
signal. Unlike the sequential algorithms, the SNR does not depend on N. Therefore 
the number of segments can be increased without suffering from noise. Like with 
the continuous sequential algorithm, the integration time can be adjusted dynam- 
ically to optimize the speed/SNR tradeoff. Although the partitioning algorithm has 
the highest SNR, the magnitude of the phase corrections decreases as the algorithm 
progresses. The required accuracy in measuring $o increases at the same pace as the 
SNR increases. 

The partitioning algorithm is very sensitive to measurement errors since, when the 
measured $o has an error, half of the segments will be programmed with the wrong 
phase. In the extreme case where the error equals n the enhancement completely 
disappears in a single iteration. A simple and effective solution to this problem is to 
keep the previous configuration of the phase modulator in memory. When an op- 
timization step causes the signal to decrease, the algorithm can revert to the saved 
configuration. 



7.6 Simultaneously optimizing multiple targets 

For some applications one wishes to focus on multiple targets at once. To achieve 
multi-target focusing, one could use the sum of the intensities in all targets as feed- 
back. In this case the feedback signal would equal (compare Eq. 7.12) 

M 

/M(0)=2]/m(0) (7.13) 

m 
M 
= ^/m^ + JmB+2v'W^COs($-$mo), (7.14) 



where the summation is over all M targets. Each of the cosine terms in Eq. (7.14) is 
weighted by the square root of the intensity contribution of the modulated segments 
ImB and the contribution of the unmodulated segments I^a- This weighting makes 
that targets that initially are more intense are optimized more effectively. While the 
intensity increases in all targets, the enhancement will not be distributed equally over 
the targets and the distribution depends on the starting conditions. 

With the sequential algorithms, it is possible to distribute the enhancement equally 
over all targets by using a different feedback signal. If we use the summed square root 
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of the intensities in the targets as feedback, we have a feedback 'amphtude' 

M 

AMW^^^/lmW (7.15) 

m 
M 

'^^\/X^+\/Z^cos($-$mo), (7.16) 



where we used that I^a "> ImB for the sequential algorithms. For these algorithms, 
ImA is approximately equal to the initial intensity in target m . The feedback signal in 
Eq. (7.16) is only proportional to the propagation amplitude ^j ImB of the modulated 
segment to the target. Therefore, the signal does not depend on the initial conditions 
and the intensity enhancement will be distributed evenly over all targets. 



7.7 Conclusion 

Three different algorithms for inverting wave diffusion were presented. The algo- 
rithms were compared experimentally, with numerical simulations and using analyt- 
ical theory. We found good agreement between experimental data, simulations and 
theory Moreover, the simulations and theory can be used to predict the performance 
in different experimental situations. 

The effectiveness of the algorithms was quantified by the enhancement. It was 
seen that the enhancement depends on the number of segments N and the relative 
persistence time Tp/Ti. For the sequential algorithms to have optimal performance, 
it is required to adjust N to match Tp. This means that these algorithms need a-priori 
knowledge of the system. The partitioning algorithm does not need this knowledge 
and always performs close to optimal. Moreover, the algorithm causes the enhance- 
ment to increase the most rapidly of the three investigated methods. All in all, this 
algorithm is a good candidate for applying inverse diffusion in instable scattering me- 
dia such as living tissue. In the future, learning algorithms (see e.g. [24, 25]) might be 
developed to further improve the performance of inverse diffusion, for instance by 
dynamically balancing the trade-off between signal to noise ratio and measurement 
speed. 

The maximum enhancement linearly depends on the number of measurements 
that can be performed before the speckle pattern decorrelates [Tp/Ti). The faster the 
measurements, the higher the enhancement. In our current system, the speed is lim- 
ited by the response time of the LCD. Fast micro mechanical phase modulators have 
a mechanical response time of about 10 [is (see e.g. [26]), which allows a 10* times 
faster operation than with our current system. At such high speeds, the response time 
of the operating system of the control computer will become a bottleneck. Since the 
computation time required for the algorithm itself is negligible, this bottleneck can be 
overcome by using a real-time operating system, or by replacing the computer with 
a programmable logic device. In perfused tissue, a typical decorrelation timescale is 
10 ms [16], which means that an enhancement of about 50 should be possible with 
currently available technology 
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Figure 7.9: Complex plane representation of the partitioning algorithm, a) Before the iteration the con- 
tributions from A and B are not exactly in phase, b) After the iteration, the contributions 
are aligned and the resulting intensity is higher. 

7.A Calculation of the performance of the partitioning 
algorithm 

In this section we calculate the development of the enhancement of the partitioning 
algorithm under ideal conditions. During one iteration of the partitioning algorithm, 
the phase modulator is randomly split into two groups [A and B), each containing 
half of the segments. The relative phase ($) of group B is cycled from to 2?!. During 
this cycle, the target intensity is given by 

m)^\E^A + EmBe'^f, (7.17) 

where E^a is the contribution of the segments in group A to the target field 



with 



/ , y ., ?mn. (7.18) 

nsA 



V Vn) 

and similar for E^b- The coefficients ^mn are initially random and distributed ac- 
cording to a normalized circular Gaussian distribution, meaning that (i^) — and 
((Ref )^) — ((Im(^)^) = 1/2. As the algorithm proceeds, the phases 0„ are adjusted 
and the distribution gradually changes to a Rayleigh distribution when a high en- 
hancement is reached. The average value (i^) increases from to y 7r/4 as all con- 
tributions are aligned to be in phase. At any moment during the optimization, {^) — 
V(??-l)/(iV-l)and(|?|2) = l. 

Figure 7.9 gives a graphical representation of a single iteration. Before the iteration, 
EmA and E^b have a different phase. Without loss of generality, we choose the phase 
oi(EmA + EmB^ to be 0. The intensity before the iteration is given by 

/before = (Re£™,,-FRe£„Bf. (7.20) 

After the iteration, <1> is set to the value that caused the highest target intensity, which 
means that E^a and E^b are now in phase. The target intensity then equals 

/after = (|£mAl+|£mBl)', (7.21) 

which is higher than or equal to before the iteration. 
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We now calculate the average intensity gained in a single iteration. We consider the 
regime where already a few iterations have been done (ry » 1). In this regime, we can 
approximate 

I £'mA I = ^/{ReEmAf + ilmE^Af (7.22) 

R^ Re£„A + ™ . (7.23) 

2ReEmA 

Using this result in Eq. (7.21) gives 

4fter = (Re^m^ + ReE^sf + {ImE^Af + (Im£'„B)^+ 

+ ^-^(Im£„Bf + — -^(ImE^^f + 

KeiimB i^ShmA 

4(Re£^^)2 4(Re£^B)' 2Re£^^Re£„B 

where the terms on the last line can be neglected. lfN:>l, then ReEmB/ReEmA '^ 1- 
The intensity gain for the iteration A7 = /after — /before is now found to be 

A/ = 2(Im£„^)2+2(Im£„B)' (7.25) 

We are primarily interested in the regime N :> r] :> I), where the algorithm picks 
up the main part of the final enhancement. In this regime, (f ) ^ 1 and the proba- 
bility distribution of f is still close to the original Gaussian distribution. Therefore, 
((Im<^)2) !%j 1/2 and it follows from Eq. (7.18) that 

A/=-(7o) (7.26) 

Therefore, we expect the intensity enhancement rj to increase with 1/2 after each 
iteration of the algorithm. With this information, we also calculate the typical phase 
adjustment that is performed in each iteration. From Fig. 7.9 it follows that the root 
mean square phase adjustment equals 



When r] approaches its maximum, all contributions are almost completely in phase 
and ((Imf )^) vanishes. In this regime, the algorithm becomes less and less effective, 
as was seen in simulations and experiments (see Fig. 7.6). Adaptive switching to one 
of the other algorithms will ensure optimal convergence. 
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Chapter 8 

Transport of light with an 
optimized wavefront 



In our experiments we construct a wavefront that maximizes the intensity at a given 
target position. In this chapter, we investigate the transport properties of light with 
such an optimized wavefront. We introduce the concept of 'active' measurements to 
describe a measurement where the wavefront is tailored to the specific realization of 
disorder in the sample. We show that this new type of experiments is very promising 
as it could be used to find signatures of absorption and possibly Anderson localiza- 
tion of light, or to test mesoscopic transport theories and selectively study extended 
modes in a sample. 

So far, we have used a simple ad hoc matrix model (Eq. (3. 1)) to describe scattering 
and focusing through a disordered medium. This model is based on the assumption 
of uncorrelated transmission coefficients (UTC). In reality, the constraint on energy 
conservation makes that the matrix elements are correlated. Such correlations can be 
described by random matrix theory (RMT). RMT predicts that it is possible to gener- 
ate an incident wavefront that is fully transmitted through an opaque object, regard- 
less of the thickness of the object. 

In Section 8. 1 we give a brief introduction to RMT and discuss how it can be applied 
to describe scattering in a disordered optical waveguide. In Section 8.2 we use RMT to 
describe the transport of light in an active measurement. We calculate experimentally 
observable quantities and compare the situations with the matched and unmatched 
incident wavefronts. Furthermore, we discuss the physical relevance of these 'actively 
measured' quantities and compare the results with predictions made with the UTC 
model. The total transmission of an optimized wavefront is found to be an observable 
that characterizes the correlations in the scattering matrix. Finally, in Section 8.3, we 
extend the RMT model so that it applies to a slab geometry and accounts for other 
experimental parameters, such as the resolution of the light modulator. We use the 
extended model to calculate how experimental limitations affect the measurements 
and we show how they can be compensated for. 



8.1 Random matrix theory 

Random matrix theory is a very general mathematical theory for statistically describ- 
ing extremely complex systems. Some excellent reviews on RMT can be found in 
Refs. (author?) 1, 2 and (author?) 3. The theory has been developed by Wigner[4] to 
describe statistics of eigenvalues of many-body systems. It was originally employed 
to describe the properties of excited states of atomic nuclei and also proved very suc- 
cessful in calculating the statistics of electron scattering in mesoscopic conductors 
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Figure 8.1 : Scattering in a waveguide geometry. The scattering matrix S couples incident waveguide 
modes £+ and E^ to outgoing waveguide modes E^, E^. In a transmission experiment, 



and quantum dots[5, 6, 7, 8]. Weaver[9] was the first to study the application of RMT 
to acoustical waves instead of quantum mechanical wave functions. In pioneering 
experiments in acoustics[9, 10] and microwaves[ll, 12] the distribution of eigenfre- 
quencies in irregular cavities was measured. RMT proved very successful in predict- 
ing the statistical properties of eigenfrequency spacing in these systems. [3] Although 
most results of RMT are for monochromatic waves and, therefore, do not describe the 
dynamics of scattering, methods for extending the matrix treatment to describe time 
dependence are being developed! 13, 14, 15, 16]. 

More recently, RMT has been used to calculate correlations of light scattering in 
disordered optical waveguides (see e.g. [2]) and random lasers[17]. Because only very 
general concepts are used, it is very easy to translate results for one type of waves to 
a different system. This also means that experimental observations in one field of 
physics are directly relevant for many other disciplines. 

The basic assumption of RMT is that the system is so complex that you might as 
well use a completely random Hamiltonian matrix to describe it. In Wigner's original 
approach[4], the 'randomness' of the matrix is invariant under unitary transforma- 
tions of Hilbert space. In terms of light scattering, this invariance means that the 
scattered light has no spatial or angular correlation with the incident light whatso- 
ever. The only other constraint that Wigner imposed on the random matrices is that 
there is a conserved quantum number that imposes a symmetry condition on the 
matrix. In the case of scattering this condition is that the energy is conserved (no 
absorption). Usually, an optical system also has time-reversal symmetry, which puts 
an extra constraint on the random matrix. Due to these constraints the elements 
of the matrix have to be statistically correlated. Typically, the number of scattering 
channels is large and the correlation between any pair of channels is negligible. Only 
when many scattering channels are involved, the effects of the correlations become 
visible. 

8.1 .1 Random matrix theory in a waveguide geometry 

We first explain how RMT is applied to a simple waveguide geometry as is shown in 
Fig. 8.1. To both sides of a strongly scattering disordered sample a waveguide sup- 
porting exactly N propagating waveguide modes (channels) is connected. Evanes- 
cent fields in the waveguide are not taken into consideration. The sample scatters 
incident light into reflected and transmitted waveguide channels. 
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Figure 8.2: Transmission experiment in a waveguide geometry. Ttie transmission matrix of ttie system 
is decomposed in three matrices: V maps the incident light to eigenmodes of the sample, 
T gives amplitude transmission coefficients for the eigenmodes. U maps the eigenmodes 
back to real space on the other side of the sample. 



Scattering by the sample is represented by the scattering matrix S that connects 
the incident field E+, E^ to the channels leaving the sample [E^ and E^). Channels 
labelled a are channels on the left side of the sample and Channels labelled b are 
channels on the right side of the sample. The scattering matrbc S is defined as 



(8.1) 



We assume that the sample does not absorb any light and, therefore, S is unitary. 
When, in addition, scattering in the sample is reciprocal, time-reversal symmetry is 
obeyed and S is both unitary and symmetric. Although aU our experimental systems 
have time-reversal symmetry, unitarity alone is sufficient for the derivations made in 
this chapter to be valid. The scattering matrbc can be split into four matrices describ- 
ing reflection and transmission 
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u\ 
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'bb 
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E7 



(8.2) 



From nov^r on, we only consider transmission from left (channels a) to right (chan- 
nels b) as described by the matrix f^^, in short t. RMT tells us that the restriction 
of energy conservation gives rise to a statistical structure in the transmission matrix. 
This structure is visible in the singular value decomposition (SVD) of the transmission 
matrix. Like any matrix, the transmission matrix t can be decomposed as 



t=UTV, 



(8.3) 



where U and V are unitary matrices and t is a diagonal matrix with real, non-negative 
elements. The physical interpretation of Eq. (8.3) is illustrated in Fig. 8.2. Matrix V 
maps the incident light from real space coordinates onto transmission eigenmodes 
of the sample. The intensity transmission coefficients for each of these eigenmodes 
are given by the diagonal elements of t^. Matrix U maps the eigenmodes back to 
spatial coordinates at the back of the sample. 

In RMT U and V are random unitary matrices with elements that have a circular 
Gaussian distribution. The statistical properties of these matrices are invariant to any 
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Figure 8.3: Probability density of tine intensity transmission eigenvalues T = t^ as derived in Ret. (au- 
thor?) 5. 



unitary transform of the Hilbert space. [4] In other words, when we replace U (or V) 
by U' — BU [V = BV), with B any unitary matrix, the ensemble properties do not 
change. It, therefore, does not matter whether a channel is defined as a waveguide 
mode, an incident plane wave or a position on the sample surface. Moreover, these 
constraints on U and V make that these matrixes are statistically independent of t. 



8.1.2 Distribution of transmission eigenvalues 

Intuitively, one might expect that all modes in a random system are roughly equiv- 
alent and the distribution of T = t^ is peaked around the total transmission of the 
sample. However, as was shown first by Dorokhov[5], the eigenvalues have a proba- 
bility density that is given by 



piT]^{T} 



1 



2TVT^^ 



(8.4) 



This probability density is plotted in Fig. 8.3. The eigenvalue distribution is bimodal: 
almost all eigenvalues are either 1 or 0. A common terminology is to speak of 'open' 
and 'closed' channels. Light in an open channel is almost fully transmitted [T i^ 1), 
whereas light in closed channels is almost completely reflected [T i^O). The predic- 
tion that all channels are either open or closed is also known as the maximal fluctua- 
tion theorem[18, 19]. 

The distribution in Eq. (8.4) is the result of the restriction of energy conservation 
that is imposed on the scattering matrix. Due to this restriction, a mechanism called 
eigenvalue repulsion spreads out the values of 1/f to a uniform distribution. [2, 18, 
20] Here, ^ is the localization length that is associated with a mode in the sample. 
The localization length relates to the transmission of a sample by T = cosh"^(L/^). 
The uniform distribution of l/^ transforms into the distribution of T as is given by 
Eq. (8.4). Originally, RMT was used to derive Eq. (8.4) for a waveguide. However, later 
it was shown that this probability density function is more generally valid and also 
applies to a slab geometry. [21] 
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Equation (8.4) has divergencies at and at 1. Tlie divergency at is not integrable. 
In fact, there is a minimum value for the transmission coefficient of a mode [2, 22], 
which is given by T = cosh"^(l/(r)). This cutoff normalizes the probability distribu- 
tion of T. The physical interpretation of the cutoff is that with (T) ^ t/L, the local- 
ization length of a mode cannot be shorter than the mean free path I. In practice, 
the cutoff for T is often ignored since most of the time one is only interested in the 
moments of T. These moments depend very little on the cutoff. [22] 

8.1.3 The effect of refractive indices 

Random matrix theory is most conveniently formulated in a waveguide geometry 
with an equal number of waveguide channels on both side of the sample. In an exper- 
iment, the refractive index in front of the sample is often different from the refractive 
index at the back. Moreover, the sample itself has an effective refractive index that 
may differ from that of its environment. The number of channels scales with the re- 
fractive index squared. Therefore, the number of channels in the sample is not the 
same as the number of channels in the connected waveguides. 

To account for the different number of scattering channels on both sides of the 
medium we assume that all channels are statistically equivalent. The number of 
channels is determined by the refractive index in each of the media. At an interface, 
scattering randomly redistributes the light over all channels. When the light enters 
a medium with a lower refractive index, some linear combinations of channels are 
completely reflected. This reflection is the random equivalent of total internal reflec- 
tion. 

In general the transmission coefficient at the surface of a disordered medium is an- 
gle dependent[23, 24, 25]. If we define channels as incident angles, different channels 
are not statistically equivalent. Therefore, we chose to define channels in real space 
coordinates. Since angular effects are equal for all channels, the channels are still 
statistically equivalent. Furthermore, we assume that the refractive index contrasts 
are low and that scattering is strong, so that neighboring channels can be considered 
statistically independent. 

We now consider three possible combinations of refractive indices and analyze 
their effect on the RMT formalism. 

Matrix description of a sample in air 

When a sample is surrounded by air, its effective refractive index is higher than that 
of the surrounding medium. This means that the sample supports more channels 
than can possibly couple to the outside world. The difference in supported channels 
is shown schematically in Fig. 8.4. Since no light can couple to these bound modes, 
we ignore them and just use the formalism in Eq. (8.3), with t a random subset of the 
actual full T matrix of the sample. 

Matrix description of a sample in a dense medium 

In the situation where the refractive index of the sample is lower than that of its en- 
vironment not all incident channels can couple into the sample. Or, more correctiy. 
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Figure 8.4: Diagrammatic representation of a sample in a medium witln a lower refractive index. The 
sample supports more channels than the connecting waveguides and can, be thought of in 
a sense being 'wider'. 




Figure 8.5: Diagrammatic representation of a sample in a medium with a higher refractive index. The 
sample supports fewer channels than the connecting waveguides and can, be thought of 
in a sense being 'narrower'. Not all incident channels can couple Into the sample. 



some linear combinations of incident channels have no support in the sample. As 
is depicted in the diagram in Fig. 8.5, the sample acts as a bottleneck. This situation 
could be modelled by using the formalism in Eq. (8.3) with a matrix t that is extended 
with extra zeros for channels that cannot exist inside the sample. 



Matrix description of a sample on a substrate 

When the sample is on a substrate with a high refractive index, it is sandwiched be- 
tween a low refractive index medium and a high refractive index medium. In this 
experimentally relevant situation, the number of incident channels is not equal to 
the number of transmitted channels and t is a rectangular M x N matrix. Here N is 
the number of incident channels and M is the number of transmitted channels. In 
the singular value decomposition of f , [/ is a unitary M x M matrix, 1/ is a unitary 
NxN matrix, and t is a rectangular diagonal MxN matrix. The number of diagonal 
elements in t equals K = min(M, AT). When the sample is illuminated through the 
low index waveguide, as is shown in Fig. 8.6, there will be no total internal reflection 
as the light travels from a low refractive index towards a high refractive index. Since 
the case with a rectangular t matrix is the most general, we will use this case in our 
derivations. 
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Figure 8.6: Diagrammatic representation of a sample on a substrate with a higher refractive index than 
that of the sample. All channels on the low index side of the sample can couple to the 
sample and, in turn, all channels in the sample can couple to the substrate. 



8.2 A new class of experimental observables 

A common way to experimentally study an optical system is to measure its response 
to a fixed incident wavefront. For a disordered system, the response will be stochastic. 
Useful results are obtained by statistical methods such as ensemble averaging and 
calculating correlation functions. 

In our experiment, however, we actively adapt the incident wavefront to each spe- 
cific realization of disorder in the sample and selectively excite an extended optical 
mode. This concept is analogous to the approach taken in the field of coherent con- 
trol, where a femtosecond pulse is shaped in the frequency domain to optimize the 
response of a chemical system[26, 27, 28, 29, 30, 31]. This pulse selectively excites a 
chemical pathway that can then be studied. For most systems described by RMT such 
control is not possible and at best ensemble averaged quantities can be observed. 

One of the potentials of our method becomes apparent from the distribution p{T] 
in Eq. (8.4). This equation predicts that there always exist channels that have 7=1. 
Using a wavefront synthesizer we can, in principle, generate an incident wavefront 
that only couples to these open channels. Such a wavefront should be fully transmit- 
ted through an otherwise opaque scattering medium. 

Open channels correspond to optical modes that extend over the whole thickness 
of the sample. [2] With our active measurement method, it is possible to selectively 
study the properties (correlations, speckle statistics, delay time statistics, behavior in 
random gain media, second harmonic generation, etcetera) of extended modes with- 
out having a contribution from the localized modes that correspond to closed chan- 
nels. The relation between localized modes and extended modes is of paramount im- 
portance in the study of Anderson localization (where the number of extended modes 
vanishes, see e.g. [2]) and in the field of random lasers[32]. 

The possibility to perform active measurements once again illustrates the funda- 
mental difference between having a light modulator in front of the sample or behind 
the sample. The two situations are shown in Fig. 8.7. In Fig. 8.7a the modulator is 
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Figure 8.7: Difference between placing a light modulator in front of the sample or behind, a) Light 
modulator in front of the sample. A sharp focus surrounded by a speckle pattern is formed 
and the total transmission in the order of 100%. b) Light modulator behind the sample. All 
transmitted light is formed into a focus without speckle. The total transmission does not 
change and typically is in the order of 10%. 



placed in front of the sample and the incident wavefront is optimized for creating a 
focus behind the sample. The optimal wavefront only couples to extended modes in 
the sample and the total transmission increases. Since only a fraction of all channels 
is open, there is no complete set of basis vectors for forming any desired field be- 
hind the sample. As a result, it is not possible to concentrate all light in the focus; the 
rest of the transmitted light will form a random speckle pattern. When a modulator is 
placed behind the sample, as is shown in Fig. 8.7b, it could in principle focus all trans- 
mitted light to a single point. However, light transport in the sample is, of course, not 
affected. 

To find a wavefront that excites only open channels, one could theoretically per- 
form a singular value decomposition of the full transmission matrix (see Eq. (8.3)). 
However, measuring the complete scattering matrix is experimentally very challeng- 
ing since it contains tens of millions of elements. Moreover, measurement noise 
would accumulate and disturb the singular value decomposition. A different, much 
easier, approach is to create a focus at any position behind the sample. Since closed 
channels are fully reflected, only open channels contribute to this focus. When we 
optimize the incident light to create a maximally intense focus, we generate a wave- 
front that only overlaps with open channels. 

In this section, we use the RMT model to calculate some experimentally observ- 
able quantities and compare the result for an unshaped wavefront (passive measure- 
ment) to the result of a wavefront that is optimized for creating a focus (active mea- 
surement). The results of the calculations that follow below are summarized in Ta- 
ble 8.1. We first derive these quantities for a waveguide geometry with a wavefront 
synthesizer that perfectly controls the phase and amplitude of all incident channels. 
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passive 


active (optimized) 


single channel 
total 







Table 8.1: Average transmission into a single channel and average total transmission with and without 
actively shaped wavefront. 

This assumption greatly simplifies the calculations. Later, in Section 8.3, we will see 
how the general results for an open system or for imperfect modulation can easily be 
found once the results for a perfect waveguide system are known. 



8.2.1 Observables in passive measurements 

We now calculate average values for the total transmission through the sample and 
for the intensity in a single transmission channel. The incident wavefront is not 
adapted to the sample, which is the normal, passive, way of measuring transmission 
coefficients. The general equation for calculating the field in outgoing channel b is 



t^b — / , ^ba t^a • 



(8.5) 



Using the singular value decomposition in Eq. (8.3), we obtain 

JV K 
Eb — / ,/ UbkT^kkVkaEa. 



(8.6) 



a k 



Because RMT is symmetric to unitary transformations of the incident field, we can 
use any arbitrary normalized field and achieve the same outcome. For simplicity, we 
illuminate only incident channel a' and have Ea — 5 a' a ■ The total transmission is 
calculated using 



M 



(8.7) 



with/fc = \Eb\^. The total transmission follows from substituting Eq. (8.6) intoEq. (8.7) 



=1 



7 UbkT^kkVka' 



After expanding the square and reordering terms, we have 

M K 
-^tot— / , / M hk^hk'^ kk^ k'k'^ka'W'a'- 



b k,k' 



Unitarity of U implies that ^^ L/^j. L/^j., — 5kk'- Therefore, Eq. (8.9) reduces to 



-S 



<kM- 



(8.8) 



(8.9) 



(8.10) 
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To find the ensemble average total transmission, we use the assumption from RMT 
that V and t are statistically independent. We can, therefore, average Eq. (8.10) over 
all realizations of disorder to find 

(rtot>= ( > jU\V,J^ ) =(t^) ( > .\Vta'\' ) =(t^) f . (8.11) 

where the unitarity of V was used in the final step. When the refractive index behind 
the sample is higher than the refractive index in front of the sample, the factor K/N 
equals one. On average, an incident field excites all eigenmodes with equal probabil- 
ity. Therefore, the average total transmission equals (t^). 

If, on the other hand the refractive index behind the sample is lower, the factor K/N 
is less than one. To understand this factor, we examine the situation where there is 
no sample, just a direct transition between the two media. When N > M, there are 
more channels in the incident field than there are in the transmitted field. The excess 
channels cannot exist in the medium behind the sample and undergo total internal 
reflection. The fraction of channels that does not experience total internal reflection 
is M/N. When a sample is present, energy is distributed randomly over all sample 
channels. Only a fraction M/N of the light can leave the sample. 

Since in the RMT formalism all transmitted channels are equivalent, the average 
intensity in a single transmitted channel is simply found by dividing Eq. (8.11) by M 

(/;,) = (t^) . (8.12) 

^ ' \ / MN 

8.2.2 Observables in active measurements 

We now calculate the same quantities as in the previous section for an active mea- 
surement where the incident wavefront is optimized to maximize transmission into 
a single channel p. We first determine this optimal incident field Ea and then use it 
to calculate the transmission into channel fi as well as the total transmission. 

Intensity in the optimized cliannel 

We now calculate the intensity in the target channel j8 after an optimization. The 
incident field that maximizes the intensity in this channel follows from the Cauchy- 
Schwartz inequality 



7 , tpaEa 



N N 

<J^\tpafJ^\Ea'f, (8.13) 

a a' 



where the two sides are equal if and only if 



'-fia' 



(8.14) 



where we denoted the optimal incident field after optimization with a tilde sign. We 
are free to choose prefactor Ao e C; we choose Ao so that Ea is normalized. Now, the 



A new class of experimental observables I 115 



optimal value of Ip follows directly from Eq. (8.13) 



ip 



= l]|f/J« 



(8.15) 



Using the singular value decomposition in Eq. (8.3), we find 

2 



JV 



K 



k 
N K 



7 UpkT^kkVka 
k 
K 

2-i2-i ^Pk^pk''^kk'^k'k' ^ca \' 



(8.16) 



(8.17) 



a k,k' 



We use the unitarity of V to get rid of the sum over a 



b-Y\^pk 

k 



■ kk- 



Our model allows us to average over U and t separately to find 

K 



/, = T^ 



M 



(8.18) 



(8.19) 



In the situation where M — K — N, the fraction of the light that is transmitted into 
channel P is equal to the total transmission before optimization (Eq. (8.11)). 

Total transmission 

We now calculate the transmitted field for the optimized incident wavefront given by 
Eq. (8.14). Using Eq. (8.6) we find 



N r K 



Unitarity of V implies that ^^ V^^ Vp^ — dkt- We now have 

K 
k 

The total transmission after optimization now follows from Eq. (8.7), 

K 



(8.20) 



(8.21) 
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b \ k J \ k' 



S' 



* -2 1 I \ u*U,„tI,„ 



(8.22) 



Using unitarity of U in the summation over b we find 



t^^oZk/ 



pk 



' kk- 



(8.23) 
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The prefactor Aq was defined to normalize the incident wavefront. Using Eq. (8.14) 
we find 

JV K 



^/3fc 



2 



^l, . (8.24) 



a k 

where the equality of Eqs. (8.15) and (8.18) was used in the last step. Finally, we have 

2 



yK 
^k 

Txot — 



Upk 



yK 



Upk' 



2 
t2 



kk 

— . (8.25) 



k'k' 



Here we see that the numerator is a weighted average over t"* and the denominator is 
a weighted average over t^ with equal weights. Averaging over U gives 

(rtot) = |:^^C4,2. (8.26) 

We see that the total transmission after optimizing the incident wavefront is equal to 
C4,2, which is defined as the ratio between the fourth and the second moment of t. 
C4,2 does not depend on the refractive indices outside the sample. Using the distri- 
bution p{T) that is predicted by RMT (Eq. (8.4)), C4,2 can be calculated 

C4,2 = -. (8.27) 

Surprisingly, C4,2 turns out to be a universal constant. A fraction of 2/3 of the power 
of an ideally optimized wavefront is transmitted, regardless of the thickness of the 
sample. 



8.2.3 Comparison with the uncorrelated model 

The RMT model is very similar to the uncorrelated transmission coefficients (UTC) 
model in the sense that both models use a matrix with statistically equivalent ran- 
dom elements. The singular value decomposition that formed the start of our RMT 
calculations can also be performed on the UTC matrix and all observables can be cal- 
culated in exactly the same way. In effect, the only difference between the two models 
is the distribution of transmission eigenvalues. When we choose (t^) to be equal in 
both models, the only observable difference is the value of C4,2. 

RMT takes into account the correlations arising from the constraint of energy con- 
servation and gives the universal result of C4,2 = 2/3. In the UTC model the matrix 
elements are uncorrelated. For such a matrix the eigenvalues are expected to have an 
exponential distribution. [33] In that case, C4,2 = 2 (t^). The power in the optimized 
channel increases to (t^), whereas the background transmission remains constant at 
(t^). The optimized transmission is twice the original transmission and, therefore, it 
is not a universal constant in this model. 

Apparently, C4,2 is an accessible observable that depends on the symmetry of the 
disordered system, which makes it very interesting to determine it experimentally. 
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Figure 8.8: Ratio of tlie fourtli moment to thie second moment of tiie transmission eigenvalues (04,2) as 
a function of sample thickness. Witfi increasing tfiickness, C4,2 rapidly converges to 2/3. [2] 

8.2.4 Random matrix theory for thin samples 

We calculated 04,2 — 2/3 using the RMT eigenvalue distribution of Eq. (8.4). This 

equation is only valid in the limit of a 'thick' sample. To investigate when a sample 

can be considered thick, we use the equation for the exact probability distribution of 

transmission eigenvalues as calculated by Beenakker[2]. The exact equation is more 

involved than its limit in the diffusion regime. Beenakker defines x = arccosh(l/T) 

and then derives 

2 
P(x,s)=-Im[/(x-iO+,5), (8.28) 



with 



4L 
and U is the solution to the implicit equation 



(8.29) 



UiC,s)^Uoi^-sUiC,s)), 



(8.30) 



with starting condition L/o(f ) — coth(f ). We used these equations to numerically eval- 
uate C4,2 and we found that €4,2 converges to 2/3 with increasing sample thickness 
(Fig. 8.8) . The convergence is very rapid: starting at C4 2 — 1 for an infinitely thin sam- 
ple, C4 2 is already at 0.7 for samples that are only one mean free path thick. There- 
fore, we conclude that samples in our experimental range can be considered thick 
and Eq. (8.4) can be used. 



8.2.5 Non-diffusive behavior 

The distribution of eigenvalues depends on the transport physics in the disordered 
medium. In the previous section, we saw that, for very thin samples, the distribution 
depends on the thickness of the sample. 

In microwave experiments, it was observed that absorption in the medium changes 
the eigenvalue distribution to an exponential distribution. [34] This result can be un- 
derstood in the sense that there is no constraint of energy conservation in such sys- 
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passive 


active 


single channel 
total 




(r,ot)=(lrP)| + (i-(lrP))(T2)f 



Table 8.2: Value of some important observables in a realistic experimental situation (also compare 
Table. (8.1)). 



terns. [33] Therefore, the elements of the transmission matrix are uncorrelated, just as 
in the UTC model. 

In the limit of extremely strong scattering {ki ss 1) interference effects dramatically 
change the transport properties of waves. [35, 36] This effect is knov^rn as Anderson 
localization. The observation of Anderson localization of light[37, 38, 39] has proven 
very difficult, especially because absorption can, in many situations, have a similar 
signature [40]. 

Both absorption and localization make that the extended modes vanish, thereby 
decreasing C4,2. However, in quasi-lD systems the effect of localization on the dis- 
tribution of eigenvalues is rather spectacular. At the localization threshold, eigen- 
values start to cluster. This effect is called eigenvalue crystallization and leads to an 
oscillating eigenvalue distribution (see e.g. [2]). The study of localization in open 3- 
dimensional systems is an ongoing field of research. Recent results[41] predict the 
distribution of decay rates of optical modes. It would be very interesting to know the 
distribution of transmission eigenmodes in such a 3-dimensional system and inves- 
tigate if, by measuring C4,2, the difference between localization and absorption can 
be established experimentally. 



8.3 Random matrix theory in an optical experimental 
situation 

In this section, we extend random matrix theory of wave scattering to describe a 
realistic optical experimental situation. First, the geometrical correlations that are 
present in a slab geometry are analyzed and then the effect of limitations of the mod- 
ulator is taken into account. We find that a realistic experimental situation can be de- 
scribed with the degree of amplitude control y that accounts for all deviations from 
the 'ideal' waveguide geometry. The value of the observables in a realistic situation 
are summarized in Table 8.2. 



8.3.1 Slab geometry 

In our experiments we have a slab geometry as is depicted in Fig 8.9. In this geom- 
etry, not all channels are equivalent. When light is focused on incident channel a, 
transmitted channels direcdy opposite a will receive more intensity than channels 
that are further away. To include this geometric effect, we recognize that the average 
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Figure 8.9: Slab geometry. Light diffuses from a point with surface coordinates ra on the front of the 
sample to a point rj, at the back of the sample. 



transmission coefficient is a function of tfie relative position of the channels. 

K 



\tba 



= /(|rfc-rJ) T^ 



MN 



(8.31) 



where all channels are defined in r space, r^ and rj, are the 2-D coordinates of, re- 
spectively, the incident and the transmitted channels at the front and back sample 
surfaces. In a waveguide geometry /(Ir^g — Tq |) = 1 . In a slab geometry / can be calcu- 
lated using diffusion theory. It is a function that starts at 1 and decreases with distance 
|r^ — To I v«th a typical decay length that is comparable to the thickness of the sample. 
Furthermore, / depends on the mean free path and the extrapolation lengths. Since 
/ decreases with distance, channels that are further away from the target contribute 
less. The effective total number of channels that contribute to the intensity in p is 
given by 



N,ff = J^f[\rp-ra\). 



(8.32) 



A/eff converges to a finite value in the limit K,M,N —>oo. In this limit, the waveguide 
model transforms into a model for an infinite slab and all physical quantities become 
independent ofN. Similarly, the effective number of transmitted channels at the back 
of the sample originating from a point source at r^' equals 



M. 



M _ 



eff = ^/(|rb-ra'|)= — Neff. 



(8.33) 



To introduce spatial dependency in the RMT formalism, we use the singular value 
decomposition (Eq. (8.3)) to find 



{\tba 



7 UbkT^kkVkc 



(8.34) 



Assuming that t is still independent of U and V, all geometrical effects are accounted 
for by imposing the following correlation between on U and V 



{\UbkVkaf)^n\rb-r„\)-^ 



MN 



(8.35) 
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In the previous sections, we (implicitly) used the assumption that the matrices U and 
V are uncorrelated only in the calculation of (lb); all other observables in Table 8T 
are not affected by this correlation. In the active measurements, the incident field 
is matched to the sample to maximize It- In a passive measurement, however, the 
incident illumination can have any arbitrary profile that only depends on the experi- 
mental configuration. The average value for lb depends on this intensity profile !„ 






When the detector is further away from the sample than the diameter of the diffuse 
spot, it is convenient to define far-field angular channels (far- field speckles). The 
number of such channels is Meff (also see Ref. (author?) 42 for a discussion on the 
number of independent channels). In that case, Eq. (8.36) reduces to 

{Ii,)^(t^) , (8.37) 

where the total incident power was assumed to be normalized. All other observables 
that were derived above are not affected by the slab geometry. It is now possible to 
supplement these results with spatial information. We can, for instance, calculate the 
intensity profile of the optimal incident field (Eq. (8.14)) 

{\E.\^) = (Al) {\tp.f) = (A^)/(|r, -r.l)(T^) ^. (8.38) 

Below, we will use this intensity profile to determine how well a generated field match- 
es the optimal field. 

8.3.2 Wavefront modulation imperfections 

So far, we have assumed that the incident wavefront is constructed with an infinitely 
large perfect wavefront generator. All practical imperfections, such as a finite size 
modulator or phase only modulation, are accounted for by a single parameter j. We 
start by decomposing the incident field into a contribution to the perfect incident 
field and an orthogonal part AEa 



Ea = rEa + Vl-lrl'A£„ , (8.39) 

where y is the inner product of the optimal field and the generated field 



JV 

r = Y.^:Ea. (8.40) 



With an ideal modulator, the optimal field can be generated perfectly and y — \- All 
experimental limitations decrease y. In the next section, we give some examples how 
to calculate y in realistic situations. Now, we first calculate how y effects the observ- 
ables Ip and Ttot. Since AEa is orthogonal to the optimal field, it does not contribute 
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to the intensity in the target focus. After optimization, the intensity in the focus is 
given by (compare Eq. (8.19)) 



{h)-m{-')l- (8.41) 

To calculate the total transmission, we use a similar approach 

(21et) = (lrOC4,2+(l-(lrP))[(T^)^-e], (8.42) 

where (t^) K/N is the transmission coefficient of the unoptimized field and e is a 
small correction term takes into account the fact that the total transmission of the 
orthogonal field is slightly lower than the total transmission of a random field since 
it does not include the contribution of transmitted channel j6 . Also, at the back of 
the sample the transmitted fields corresponding to /S.E and Ea are not completely 
orthogonal. Both effects are of order 1/A/eff and partially cancel. Therefore it is safe to 
neglect e. 

When the incident field is not shaped to fit the scattering matrix of the sample, 
there will still be some accidental overlap with the optimal field. For such an incident 
field, (I7P) is of order l/Ngg and will be neglected. 

8.3.3 Examples of realistic experimental situations 

We now calculate (li'p) for various experimental conditions. Eqs. (8.41) and (8.42) 
can then be used to obtain the intensity in the focus, as well as the total transmission 
after optimization. 

Example 1 : Finite size modulator 

Suppose the light modulator illuminates a finite area that contains Ns independent 
scattering channels. The modulator perfectly modulates amplitude and phase of all 
these scattering channels. Then 

N, JVs 



"^ 



ElEa^AoAfJ^\tpa\\ (8.43) 



where Ao normalizes the optimal field (see Eq. (8.24)) and Af normalizes the actually 
generated field. Using Eq. (8.31) to introduce position dependence, we find 

^'^ r.fi\rp-ra\) ^eff ■ 

In the limit where the light modulator is very large j ^> I. In the limit where the 
modulator is much smaller than the diffuse spot described by /, the degree of control 

is(|7-|2)=JV,/iVeff. 
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Example 2: Limited resolution 

We first calculate (I7P) for the situation where the light modulator has a limited spa- 
tial resolution. We assume that each of the Ns segments of the phase modulator il- 
luminates a spatial cluster of N/Ns scattering channels. We group the channels by 
segment and use Eq. (8.40) to find 



r 



=11 



E*E^ 



(8.45) 



where {as} is the collection of channels that is illuminated by segment s. £5 is the 
field that is generated by segment s of the modulator. The incident field is optimized 
for a whole cluster at a time. After an optimization we have 



= A X 



As = 



N 






X Ea 



-1/2 



(8.46) 



(8.47) 



where the prefactor As normalizes the total power of the incident field. Substituting 
Eq. (8.46) in Eq. (8.45) we find 



-A, 2] 



S 



^r^^n 



X Ea 



aG{a^\ 



Because elements of £„ are independent and Ea is normalized, we have 

1 



I ^a^a 



'■5 a', 



N 



which leads to 



m-^- 



(8.48) 



(8.49) 



(8.50) 



(8.51) 



Example 3: Phase only modulation 

We consider the case of a phase only modulator that perfectly controls the phase of 
all N incident waveguide channels. The amplitude of the channels is fixed at 1/vTV. 
For this experimental configuration, we find 



JV . N 

L^ a a JJj/—l^ °" 
a ^ a 

1 ^ ^ ^ 1 

- X \Ea\\Ea'\+-Y:Ea"?. 



(8.52) 



(8.53) 



a,a'^a 
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To average, we use the complex Gaussian probability function of Ea to arrive at 

Example 4: Realistic experimental conditions 

In this example, we calculate the degree of control under realistic experimental con- 
ditions. We analyze the geometry that was used for measuring the data in Fig. 3.4 
where the enhancement as a function of the number of modulator segments is plot- 
ted. In this geometry, no lenses were used behind the sample; the camera was ap- 
proximately 30 cm away from the sample surface. The incident light is focused on the 
sample surface and illuminates a spot that is much smaller than a millimeter. Analy- 
sis of the measured tba coefficients showed that, on average, all segments contribute 
approximately equally to the amplitude of the target (variations are within 15%). Fur- 
thermore, limitations of the phase modulator made that the generated field had a 
bias. On average, 14% of the incident intensity had a constant phase and could not 
be shaped. 

Since all segments contribute equally, we expect the enhancement to be given by 
a combination of Example 2 and Example 3. Also taking into account that 14% of the 
incident light could not be shaped, we find 

(v 2) =0.86^ + , (8.55) 

which is an upper limit for the case where the sample is absolutely stable and there 
is no influence of measurement noise. We can now calculate the expected maximum 
intensity enhancement r] 

r]=j^- (8.56) 

W 

Using Table 8.2 with Eq. (8.55) yields 

r]^{\rf)M,g- (8.57) 

= 0.68(W,-1) + 1. (8.58) 

This result is the same as was found with the UTC model in Chapter 3. When all seg- 
ments contribute equally, the UTC model gives the correct results even if the sample 
is not in a waveguide geometry. The measured enhancement (17 ^ 0.57 Ns , see Fig. 3.4) 
was slightly lower than the theoretical value because of measurement noise and due 
to the fact that segments in the center contribute slightly more than segments close 
to the edge of the phase modulator. This minute effect could be included in / in 
Equation (8.44), to calculate a value for (Ifp) that is a few percent lower. In the exper- 
iments described in Chapter 3 the background intensity did not increase noticeably. 
This result is explained by the fact that Ns <^ Ngg and, therefore, (Ij'p) <i I. In this 
regime, the total transmission is hardly affected by an optimization (see Eq. (8.42)). 
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8.4 Conclusion 

We developed a model for describing the transport of a waves that are actively shaped 
to focus through a disordered sample. The new model is based on RMT and, there- 
fore, correctly takes into account correlations between scattering channels. 

The concept of 'active' measurements was introduced. These measurements are 
based on wavefront shaping. One of the most interesting observables in an active 
measurement is the total transmission of an optimized wavefront. The fraction of 
transmitted power for such a wavefront turns out to be 2/3, which is universally true 
regardless of the thickness of the sample. Physical processes such as absorption and 
Anderson localization change this value, which makes measuring it extremely inter- 
esting. 

All results were first calculated in an ideal waveguide geometry and then translated 
to a realistic experimental slab geometry. This two-step calculation process proved 
to be a very powerful and versatile technique. 

Apart from the observables that were calculated in this chapter, there is a multitude 
of interesting quantities that can be measured with our active measurement strategy. 
Since most of the theoretical work has been done on ensemble averages rather than 
on the subset of extended modes, many questions are still unanswered. 
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Chapter 9 

Observation of open transport 
channels in disordered optical 
systems 



When a beam of light impinges on a strongly scattering object, such as a layer of white 
paint, tiny particles randomly scatter light in all directions. After a few scattering 
events, light is completely diffuse. The transport of multiply scattered light is, in gen- 
eral, described very well by a diffusion equation. However, the diffusion equation 
does not take into account interference effects. The inclusion of interference in the 
description of disordered wave propagation has led to the discovery of fundamental 
effects like enhanced backscattering[l, 2], universal conductance fluctuations [3, 4] 
and Anderson localization of light[5, 6, 7]. 

Here we demonstrate experimentally a new interference effect that was first pre- 
dicted by Dorokhov[8]. Like the interference effects mentioned above, this effect was 
first derived for the propagation of electron wave functions in a disordered waveg- 
uide (a metal wire) and then translated to optics. [9, 10] A disordered waveguide can 
be described as a collection of orthogonal transport channels that each carry a frac- 
tion of the incident wave. Using random matrix theory (RMT) , Dorokhov showed that 
due to interference these channels are either open (transmission close to 1) or closed 
(transmission close to 0). 

The existence of open transport channels means that it is, in principle, possible to 
construct a wave that is fully transmitted through a normally opaque mesoscopic ob- 
ject, regardless of the thickness of that object. In a mesoscopic electronic conductor 
such control over the incident wave is not possible. For light, however, we can shape 
the wavefront of the incident wave with a high accuracy. In this chapter we show that 
the shaped wave selectively couples to open transport channels, causing the trans- 
mission through an opaque, strongly scattering sample to increase dramatically. 

In Section 9.1 we give the theoretical relation for the total transmission of a shaped 
wavefront. Then, in Section 9.2, the wavefront optimization setup is described. The 
experimental results are presented in Section 9.3, followed by a more detailed de- 
scription of the data analysis procedure and an estimation of the experimental error 
margins in Section 9.4. 



9.1 Expected total transmission 

The average transport of light is described accurately by a diffusion equation. On 
average, the total transmission coefficient is given by [I + Zei)l{L + Ze\ + 2^2), with i 
the mean free path for light, L the thickness of the scattering sample and Ze\ and Ze2 
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extrapolation lengths that account for reflection at the sample interface[ll, 12]. 

To increase the total transmission above the value predicted by diffusion theory, 
we use a feedback algorithm that shapes the incident wavefront. The algorithm max- 
imizes the intensity in a target focus behind the sample. Since only open transport 
channels contribute to the intensity in the target, the ideally shaped wavefront se- 
lectively couples to these open channels. In the extreme case where all open chan- 
nels have a transmission coefficient of exactly 1 (the so called maximum fluctuations 
theorem[13, 14]), light with an optimal wavefront is completely transmitted. 

Random matrbc theory[8] predicts a distribution that is close to the extreme situa- 
tion of maximal fluctuations (see Eq. (8.4)). We express the amount of fluctuations in 
the distribution as C4,2, which is defined in Eq. (8.26) and denotes the transmission 
coefficient of an ideally shaped wavefront. For Dorokhov's distribution C4,2 = 2/3. 
This universal value is independent of the thickness or the mean free path of the sam- 
ple. 

To couple to open channels, a high degree of control over the incident wavefront 
is essential. When only a small subset of the scattering channels is controlled, the 
transmission matrix elements are effectively uncorrelated and no open channels can 
be found. Therefore, the experiment is designed to approximate the ideal wavefront 
as closely as possible. Still, in any experiment there will be a difference between the 
optimal wavefront and the field that is actually generated. This is especially true in an 
open system where the optimal wavefront has infinite support. The total transmis- 
sion of an incompletely optimized wavefront equals (see Eq. (8.42)) 

f,o,= ^ = lrl'C4,2+(i-lrl')rtot, o.i) 

where Ptot is the totally transmitted power of the incompletely optimized wavefront. 
Pin is the incident power and Ttot = Ptot/Pm, with Ptot the totally transmitted power of 
an non-optimized wavefront. It^P is the degree of control as defined in Eq. (8.40). In 
the case that we have full control over the incident wavefront Ij'P ^ 1 and Ttot —> C^^z. 
The degree of control can be determined experimentally by measuring the power 
in the target focus Pp after optimization. From Eq. (8.41) we can make \x\^ explicit 

lrP = /^, (9.2) 

where ni — \ and ^2 — 1-52 are the refractive indices of the medium in front of the 
sample (air) and behind the sample (glass substrate) respectively. Since \y\^ can be 
determined experimentally, Eq. (9.1) can be solved to find C4,2. 



9.2 Experiment 

The experimental apparatus is designed to have optimal control over the incident 
wavefront. The setup is shown in Fig. 9.1. An expanded beam from a 632. 8nm HeNe 
laser is rotated to a 45° linear polarization by a half waveplate and impinges on a 
polarizing beam splitter cube. Both horizontal and vertical polarization are modu- 
lated using the macropixel modulation method that was described in Section 2.1.4. 
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Figure 9.1 : Experimental setup for finding open transport cliannels. HeNe, expanded 632.8 nm HeNe 
laser; A/2, half waveplate; PBS, polarizing beam splitter cube; SLM, spatial light modulator. 
A, iris diaphragm; 63x, microscope objective; 60x, oil immersion microscope objective; 
S, sample; P, polarizer; ND, neutral density filter; Li, L2, L3, lenses with focal length of 
respectively 250 mm, 200 mm and 600 mm. ND and S are translated by computer controlled 
stages. 



The synthesized wavefront is focused onto the sample with a high numerical aper- 
ture (NA) objective^ so that almost all angular channels are controlled. A high NA 
oil-immersion objective^ collects almost aU transmitted light. The immersion oil 
matches the refractive index of the substrate n — 1.52 and, therefore, eliminates re- 
flections at the back surface of the substrate. The transmitted light is split into hori- 
zontal and vertical polarizations by a beam splitter cube. A second polarizer improves 
the cube's limited extinction ratio for reflected light. A camera measures the power of 
the horizontally polarized light in an area with the size of a single speckle to provide 
feedback for the optimization algorithm. We use the stepwise sequential algorithm 
(see Section 7.1.1) to optimize the wavefront. A calibrated neutral density filter with 
a transmission of 1.4-10"^ is automatically placed in front of the camera to be able to 
measure the high intensity in the target focus after optimization. The magnification 
of the detection system is 225 x. At this magnification, individual speckles are well 
separated on the camera. A second camera measures the transmission in the verti- 
cal polarization channel. This camera is not used during the optimization procedure 
itself. 

The sample is a layer of spray-painted ZnO particles (see Section 2.4) in an air ma- 
trix. The particles have an average diameter of 200 nm, which makes them strongly 
scattering for visible light. The mean free path was determined by measuring the to- 
tal transmission and equals 0.85 ±0.15 p,m at a wavelength of 632.8 nm. Two samples 
were used: one with a thickness of 5.7 ± 0.5 jim, and one that is approximately twice 
as thick (11.3 ± 0.5 nm). The sample was mounted on a motorized translation stage 
that translates the sample perpendicular to the axis of the objectives between exper- 
iments. This way, a different random configuration of scatterers could be used for 
each experiment. The translation stage also provides a reliable way of making refer- 
ence measurements (see Section 9.3). 



'Zeiss 440068-0000-000 Achroplan 63x/NA = 0.95 without cover glass correction. 
^Nikon MRD01691 CFI Plan Apociiromat TIRF 60x/NA= 1.49 
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Figure 9.2: Logarithmic piot of tlie intensity transmission in the horizontal polarization channel. The 
intensity is integrated per column, a) 11.3|im-thick layer of ZnO pigment. Solid curve, 
transmission of non-optimized wavefront; Dashed curve, transmission of optimized wave- 
front. Measurements with and without neutral density filter were combined to cover the full 
dynamic range, b) Same for a 5.7|im-thick layer of ZnO pigment. 



9.3 Results 

Optimizing the incident wavefront caused tlie intensity in the target to increase dra- 
matically. The increase in the target intensity was a factor of 746 ±28 for the 11.3 \im- 
thick sample and a factor of 671 ± 14 for the 5.7^xm-thick sample in the particular 
cases that are shown in Fig. 9.2. In addition, the intensity also increased outside the 
target area, even though the algorithm only used the intensity in the target as feed- 
back. This is an important observation since it indicates that the transmission matrix 
is correlated. By optimizing transmission to a single target, the incident intensity is 
redistributed from closed channels to open transport channels. If the elements of 
the transmission matrix had been uncorrelated, the background intensity would not 
have changed (also see Section 8.2.3). 

The increase in the background intensity is limited to an area of a few micrometers 
around the target. The diameter of this area seems to be roughly proportional to the 
thickness of the sample. This makes sense as, in a slab geometry, propagation paths 
that are further apart than the thickness of the sample are essentially uncorrelated. 
The intensity increase of the background is maximal just around the optimized spot. 
There, the intensity increased by a factor of 4.1 ±0.7 and 2.3± 0.3 for the thick and the 
thin sample respectively^. The total transmission in the background (not including 
the optimized spot) increasedafactor of 1.32 ±0.05 and 1.26 ±0.08, respectively. 

The algorithm optimizes the transmission through the whole optical system that is 
between the light modulator and the camera. To make sure that the extra intensity 
is indeed the result of an increased transmissivity of the sample, we need to exclude 
the trivial effect of an increased transmission trough e.g. the microscope objective. 



^In Fig. 9.2 the intensity increase appears to be lower since the figure shows the intensity integrated over 
all columns. 
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Translation (i^m) 

Figure 9.3: Typical total transmission after optimization as a function of sample displacement. Solid 
curve, horizontal polarization. Dashed curve, vertical polarization. Curves are normalized 
to the reference transmission (average intensity between 0.6 |im and 3.9 |im displacement). 
Fluctuations of the intensity are the result of speckle. 



Therefore, we used exactly the same wavefront for measuring the non-optimized ref- 
erence transmission as we used to measure the optimized transmission. The only 
difference between the two measurements is a slight translation of the sample. At 
the new position, the scattering behavior of the sample is different and the optimized 
wavefront has no special significance anymore. 

The effect of translating the 11.3 jxm-thick sample is shown in Fig. 9.3. Initially, the 
total transmission in both of the polarization channels has increased significantly. 
The intensity on one of the cameras is slightly higher since it includes the contribu- 
tion of the optimized spot. Since the pixels in the target area saturate, in this plot 
the difference is relatively small. After translating the sample over 600 nm, the effect 
of optimizing the wavefront is completely lost. The transmission is averaged over 34 
images taken at displacements increasing from 0.6 fim to 3.9 [im to obtain the refer- 
ence transmitted power for the non-optimized situation (Ptot)- We observe a slight 
variation in the total transmission when the sample is translated. The expected vari- 
ance of fluctuations in the total transmission is l/(2g) (see e.g. Ref. (author?) 15), 
where g = NTtot, with N the number of participating scattering channels. From the 
size of the diffuse spot we estimate N i^ 4.5 • 10^. With a measured transmission of 
Ttot — 0.24 (see Section 9.4), we expect fluctuations with a standard deviation of 2%, 
which is consistent with the data. 

We performed automated sequences of optimizations. Between each measure- 
ment, the sample was translated to find a new realization of disorder. For each of 
the measurements, we determined the degree of control using Eq. (9.2) . In Fig. 9.4 we 
plotted the increase in total transmission versus the degree of control for a 11.3p,m- 
thick sample. The total transmission increases linearly vnth Ij'p. The maximum in- 
crease in transmission that was achieved is 40%. We conclude that, with increasing 
degree of control, a larger fraction of the incident intensity couples to open channels, 
causing the total transmission to increase. 

From a linear regression we find that the curve has a slope of 1.7. With an original 
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Figure 9.4: Normalized total transmission as a function of the degree of control lyp for a single se- 
quence of measurements. Circles, measured transmission after constructing the incident 
wavefront; solid line, linear regression to the measurements; dotted line, diffusion theory; 
dashed line, random matrix theory, €4,2 =2/3 (there are no adjustable parameters). 



average transmission coefficient of Ttot — 0.24 we expect a total transmission of 0.65 
when \j\^ — 1. Tliis result is very close to the expected value of 2/3. 

For this particular sequence, we reached a maximum control fraction of |)-p = 0.21. 
The algorithm used 3816 control segments and phase only modulation. In this situ- 
ation Ij-p can theoretically not be higher than (3816/Ar)(7r/4) iv 0.67. Even in a per- 
fecdy stable system, lyf will be lower due to geometrical effects (see Section 8.3.2), so 
Ij'p — 0.21 is in line with expectations. The variations in Ij-p are probably the result of 
drifting environmental conditions that affect the stability of the setup. Although this 
drift is undesired, it has the advantage of giving a vnde range of |)'p to investigate. 

We now combine the results from different measurements into a universal plot. To 
do this, we subtracted the contribution of the uncontrolled fraction of the incident 
wavefront to find the controlled part of the transmission 



Tr 



■(i-lrl')rto 



(9.3) 



FromEq. (9.1) we expect that T^ = 17^04,2. As a result, all measurements should follow 
a universal curve with a slope of €4,2. In Fig. 9.5 we plotted the controlled transmis- 
sion for seven measurement sequences performed on two different samples. All data 
points collapse to a single line, regardless of the thickness of the sample. This is a clear 
indication that the transmission of the ideal wavefront (04,2) is a universal parameter. 
For the 11.3p,m-thick sample the linear regression gives a slope 0.68 ±0.01 and 
for the 5.7txm-thick sample we find a slope of 0.68 ± 0.02, where the uncertainties 
represent the 0.95% confidence interval of the regression. Combining the data from 
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Figure 9.5: Total transmission of the controlled fraction of the wavefront. Data from seven measure- 
ment sequences on two different samples collapses to a universal curve. There are no free 
parameters in the data or the theory. Circles, results for 1 1.3 |im-thick sample; squares, re- 
sults for 5.7 nm-thick sample; solid line, linear regression with slope of 0.68. The theoretical 
curve for C4,2 = 2/3 is not shown as it overlaps with the fit. 



both samples, we find C4,2 = 0.68 ± 0.07 where the uncertainty is determined by the 
worst case systematical errors due to the extrapolation of the intensity profile (see 
Section 9.4). 



9.4 Details of the data analysis 

In this section we give experimental details about the data analysis. In Section 9.4.1 
we explain how the totally transmitted power was measured using a CCD camera. 
This measurement method introduces a small but significant error in the determina- 
tion of C4,2. Other possible sources of systematical measurement error are discussed 
in Section 9.4.2. In Section 9.4.3 we explain how the incident intensity was measured 
accurately. 



9.4.1 Diffuse transmission measurement with a camera 

We measured the totally transmitted power with the CCD cameras behind the sample. 
The cameras collect most of the transmitted light, but the diffuse spot has such a wide 
profile that its tails are not imaged onto the CCD chip. Therefore, the camera image 
is extrapolated to obtain the total intensity. 

It turned out to be hard to reliably fit a theoretical curve to the complete intensity 
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Figure 9.6: Logarithmic plot of the transmitted intensity profile for a reference measurement. The inten- 
sity is integrated per column. Dots, measurements; Dashed line, exponential extrapolation; 
Dotted line, steepest slope that is consistent with diffusion theory 



profile, as the profile depends on the incident wavefront and many other parameters. 
According to diffusion theory, the transmission profile asymptotically decays expo- 
nentially with a decay length of (L + Zgi + Zez)!'/^ (see Section 1.3.4). It, therefore, 
makes sense to extrapolate the camera image exponentially. 

A typical transmission profile of a reference measurement is shown in Fig. 9.6. It 
can be seen that the exponential extrapolation gives a slight overestimation of the 
intensity since the measured profile does not fully reach into the exponential asymp- 
totic limit. Alternatively, we extrapolated the profile with the steepest exponential de- 
cay that is in agreement with diffusion theory. For this, we took the lower limit for the 
sample thickness (L — 10.8 [xm and L — 5.2 p,m), the minimal mean free path [l — 0.7) 
and the minimal extrapolation lengths (z^i — Ze2 — 2/3£). This alternative extrapo- 
lation method (dotted lines in Fig. 9.6) gives a lower limit for the total intensity. By 
averaging the results of the two methods we find a total intensity of (2.92 ± 0.08) • 10^ 
counts for the profile in Fig. 9.6, where the error margin is the difference between the 
results of the two extrapolation methods. A fraction of 82% of the total intensity fell 
on the CCD. 

The typical uncertainty in this experiment is about 4%. The profiles for the ref- 
erence intensity Ptot and the optimized intensity Ptot are almost equal at the edge of 
the image (see Fig. 9.2). Therefore, it is likely that the error caused by the extrapo- 
lation is either positive for both profiles, or negative for both profiles. In this case, 
the extrapolation only causes an uncertainty of ±2% in C4 2. By calculating 04,2 for all 
combinations of error margins a worst case error margin of ±9% is found. 



9.4.2 Possible causes of systematical error 

It is possible that small systematical errors are introduced by the imaging system (mi- 
croscope objective, lenses and camera) that record the transmitted light. We mea- 
sured the intensity distribution coming from a sample that is illuminated with a wide 
(> 500|jLm) spot. In this situation, the intensity distribution of the transmitted field 
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should be flat (except for speckle). We used this flat-field reference image to compen- 
sate for the position dependent sensitivity of the imaging system. 

The angular dependency of the imaging system was not calibrated. However, the 
width of the optimized spot was equal to the diffraction limit of A/(2n2) to within an 
experimental uncertainty of 6%. Therefore, we expect the effect of angular depen- 
dency to be below 6%. 

9.4.3 Measurement of the incident power 

With the cameras we measure the total power that is transmitted through a sample 
before and after an optimization. We also measure Pm, which is a reference for the 
incident power. This reference value is equal to the total power that is incident on 
the sample times the detection efficiency of the transmitted light. When Pm is known, 
we can calculate absolute values for the transmission coefficients. In measuring Pm, 
great care was taken not to introduce systematical errors. A complicating factor is 
that the optimized wavefront is different for every measurement, which causes vari- 
ations in the power that reaches the sample. To compensate for these variations, we 
determined P,n using 

P. = P..|^ (9.4) 

where Pi is the detected transmitted power with an empty substrate in place of the 
sample. The substrate was illuminated with an unshaped wavefront: all pixels of the 
modulator were set to the same phase, so that the incident light forms a focused wave. 
Pz was measured using the same unshaped wavefront, with the sample in place. The 
sample was at exactly the same position where the optimization was performed. Ptot 
is the reference transmission of the shaped wavefront. It was measured by translating 
the sample after each optimization (see Fig. 9.3). The factor 1/(1 — i?) compensates for 
the Fresnel reflection at the front surface of the substrate. For a uniformly illuminated 
microscope objective with NA — 0.95, the average Fresnel reflection coefficient is R — 
0.11. We determined the transmission coefficient of the unshaped wavefront using 



9.5 Conclusion 

We have found experimental evidence of the existence of open transport channels 
for light in opaque, strongly scattering materials. A wavefront shaping algorithm was 
used to focus light through opaque objects. We observed that the transmitted in- 
tensity in the target focus increased by almost three orders of magnitude. Simulta- 
neously, the intensity outside the focus increased, even though the algorithm only 
optimized the intensity in the focus. This result cannot be explained by diffusion 
theory or by an uncorrelated matrix model. We conclude that the elements of the 
transmission matrix are correlated as predicted by random matrix theory. Optimiz- 
ing the transmission at a single point creates a wavefront that selectively couples to 
open transport channels in the medium. 

The total transmission has a linear relation to the intensity in the target focus. This 
relation confirms the theoretical results from Section 8.3.2. The total transmission 
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through a disordered sample was increased by a factor of up to 1.4. 

We used the intensity in the target focus to derive the degree of control over the 
incident v^ravefront. By considering only the transmission of the controlled fraction of 
the v^ravefront, we were able to compare the results for different samples and different 
realizations of disorder. All results showed a universal behavior that is in excellent 
quantitative agreement with the distribution of channel transmission coefficients as 
is predicted by random matrix theory. 
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Chapter 10 

Summary and outlook 



In this thesis, we explored the use of wavefront shaping to steer light through strongly 
scattering materials. We found that scattering does not irreversibly scramble the inci- 
dent wave. By shaping the incident wavefront, an opaque object can form any desired 
wavefront and focus light as sharply as an aberration free lens. We used feedback 
from a target behind, or in, an opaque object to shape the incident wave. This way, 
light was focused through, or inside, an opaque object for the first time ever. 

Wavefront shaping provides a new way to perform experiments on strongly scatter- 
ing samples. By studying the total intensity transmission of a shaped wavefront, we 
found experimental evidence of open transport channels for light. This is an example 
of how wavefront shaping can be used to study fundamental properties of scattered 
waves. 

Our results lay both an experimental and a theoretical basis for wavefront shaping 
experiments. The concepts that were developed are generally applicable to any linear 
complex system. Therefore, our methods can be directly translated to other disor- 
dered systems such as scattering of microwaves, sound, seismic waves, electrons or 
neutrons. Even though we have only taken some first steps in this field, the outline of 
future applications and exciting research subjects becomes visible. We now discuss 
some of the possible applications of our work. 

Opaque object can be made to focus light with an excellent resolution. Therefore, 
thin layers of paint might be used as cheap disposable lenses. Scattering can also 
be used to improve the resolution and, thereby, the data transfer rate of a free space 
communication system. Instead of creating a focus, multiply scattered light can also 
be made to form a beam. This could be useful when performing laser spectroscopy 
of opaque materials. The scattered coherent light could then be collected with a high 
efficiency, thus greatly improving the sensitivity. Even better, incoherent background 
light does not form into a beam and, therefore, would not be detected. The wavefront 
that focuses light through an opaque object is unique for each object. This property 
could be used in security applications. When, for example, a small photodiode is 
covered with a layer of white paint, it is only possible to focus light on it when the 
matching wavefront is known. 

Light can also be focused inside an opaque object. The applications of this method 
are endless. One could, for example, focus light on a single fluorescently labelled cell 
below the human skin or inside bone. This application would allow one to study, 
track or destroy individual (malignant) cells. To make this futuristic application pos- 
sible, the wavefront needs to be constructed faster than the cells in the tissue move. 
We demonstrated that wavefront shaping in a living flower petal is possible. For 
wavefront shaping in human tissue, especially in perfused tissue, the speed of the 
wavefront shaping process needs to improve dramatically. This speed improvement 
is possible with the use of extremely fast light modulators that have recently been 
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developed. 

Our method can also be used to focus light inside non-random complex systems. 
This way, one can focus light on a quantum dot in a three-dimensional photonic crys- 
tal cavity, or excite a desired mode in a complex photonic metamaterial or plasmonic 
device. All in all, v^ravefront shaping can be used to guide light (or any other wave) to 
any desired position in an arbitrary complex environment. 



Nederlandse samenvatting 



Dit proefschrift gaat over de voortplanting van licht in ondoorzichtige witte dingen 
zoals velletjes papier, laagjes verf, eierschillen, bloemblaadjes en tanden. Witte ma- 
terialen verstoren de voortplanting van licht: als je bijvoorbeeld met een zaklantaarn 
op een velletje papier schijnt, gaat het licht niet als een bundel door het papier been, 
maar in plaats daarvan gloeit het hele papier op met een vaag schijnsel. Het licht is 
diffuus geworden: het straalt gelijkmatig naar alle richtingen uit en er is niet meer te 
zien v^^aar het licht oorspronkelijk vandaan kwam. Het is daarom ook onmogelijk om 
door het materiaal been te kijken. 

Tot nu toe was het niet mogelijk om te beinvloeden hoe licht zich door een on- 
doorzicbtig materiaal voortplant: het licht v^rerd altijd diffuus. Wij hebben een tecb- 
niek bedacht v^raarmee we het licht wel kunnen sturen. Met deze techniek kun je 
licht door een ondoorzichtige laag been projecteren of zelfs concentreren op een 
punt diep binnenin een ondoorzichtig object. Hiermee wordt bet misschien ooit 
mogelijk om licht door buid of hot been op een enkele gemarkeerde kankercel te fo- 
cusseren, of om door de ingewikkelde structuur van een microchip been een enkele 
transistor te belichten. Ook bij andere toepassingen, zoals telecommunicatie, optis- 
che experimenten aan microstructuren en spectroscopic aan poeders kan onze tech- 
niek uitkomst bieden. Uiteindelijk zou dit ook kunnen leiden tot een verbetering van 
bestaande medische beeldvormende technieken. 



Interferentie en verstrooiing 

Licht is een golf Zoals aUe golven kunnen lichtgolven met elkaar interfereren. In- 
terferentie is niets anders dan dat twee golven bij elkaar opgeteld worden. Als twee 
golven tegelijkertijd aankomen zijn ze in fase, zoals te zien is in Figuur la. De pick 
van de ene golf telt dan op bij de pick van de andere golf en het resultaat is een golf 
met grote pieken en dalen (dit beet constructieve interferentie). Het kan ook zijn dat 





Figuur 1 : Interferentie. a) De golven komen tegelijk aan. Ze zijn in fase en versterken elkaar. b) De 
goiven zijn uit fase en doven elkaar uit. 
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Figuur 2: Principe van het experiment, a) een vlakke golf wordt gefocusseerd op een ondoorzichtig 
wanordelijk materiaal. Door verstrooiing ten gevolge van microscopisclie oneffenlieden legt 
het iiclnt een zeer complexe route at. Na honderden keren verstrooien bereikt een klein deel 
van het licht de achterkant van het materiaal en vormt een diffuse viek. b) Door de fase van 
het binnenkomende licht op precies de goede manier aan te passen is het mogelijk om de 
bijdragen van alle vlakjes constructief te laten interfereren op een punt achter het medium. 
Het medium werkt in dit geval als een wanordelijke lens en focusseert het licht. 



twee golven uit fase zijn. De ene golf komt dan wat later aan dan de andere zodat de 
pieken van die golf samenvallen met de dalen van de andere golf; het resultaat is dan 
dat de twee golven elkaar opheffen. Deze situatie is afgebeeld in Figuur lb. 

In een wit materiaal zitten allemaal kleine deeltjes waar het licht op botst. Bij iedere 
botsing verandert de voortplantingsrichting van het licht. Het licht 'stuitert' door het 
materiaal, een beetje zoals het balletje in een flipperkast. Dit proces heet verstrooi- 
ing van licht. Al het verstrooide licht interfereert met elkaar, zodat er een complex 
wanordelijk interferentiepatroon ontstaat. Dit patroon staat bekend als laserspikkel 
(Engels: speckle). Een voorbeeld van een spikkelpatroon is te zien in Figuur 3.3a in 
Hoofdstuk 3. Het hangt af van de fase van alle interfererende golven of ergens een 
heldere of juist een donkere plek ontstaat. 



Licht focusseren met een eierschil 

We hebben een opstelling ontworpen die een bundel laserlicht in duizenden vlak- 
jes verdeelt en dan de fase van elk van deze vlakjes afzonderlijk kan aanpassen. Het 
licht van elk vlakje zal op een andere manier verstrooid worden (zie Figuur 2a) . Door 
de fase van de vlakjes op de juiste manier in te stellen kunnen we nu zorgen dat het 
verstrooide licht op een bepaalde plek constructief interfeert (Figuur 2b). De licht- 
intensiteit op die plek neemt daardoor enorm toe, soms wel met een factor duizend. 
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Het resultaat is dat een ondoorzichtig object, zoals een eierschil, licht focusseert tot 
een scherp punt. Eigenlijk werkt de eierschil dus als een 'ondoorzichtige lens'. 

Er zijn ongelooflijk veel mogelijke combinaties van fases te maken. Om de juiste 
combinatie te vinden gebruiken we een computer die de hoeveelheid licht in het 
gewenste focus meet en de fases aanpast. In het begin is de lichtintensiteit heel laag, 
maar door de fase van een vlakje aan te passen kunnen we de hoeveelheid licht laten 
toenemen. De computer herhaalt dit proces voor alle vlakjes. Uiteindelijk hebben we 
dan een heel ingewikkeld gevormde bundel gemaakt die precies compenseert voor 
alle botsingen aan de deeltjes in het materiaal. In hoofdstuk 3 laten we zien dat deze 
methode de optimale combinatie van fases vindt. Omdat deze combinatie afhangt 
van de positie van alle deeltjes in het materiaal, werkt onze methode nog niet voor 
materialen waarin de deeltjes bewegen, zoals melk of doorbloed weefsel. Met recent 
ontwikkelde technologie is het echter mogelijk om het proces zo te versnellen dat we 
in de toekomst 'real-time' kunnen compenseren voor de bewegende deeltjes. 

Belangrijkste resultaten 

Het belangrijkste resultaat van dit proefschrift is dat we in staat zijn om licht door 
een ondoorzichtig materiaal been te focusseren. Deze doorbraak heeft nationaal en 
internationaal veel aandacht gekregen (zie bijvoorbeeld 'Witte verf of eierschaal als 
hrandglas voor bundel laserlichf in het NRC Handelsblad van 14 juli 2007). Het is niet 
alleen gelukt om verstrooid licht in een scherp focus samen te laten komen, we be- 
grijpen de resultaten ook kwantitatief. We hebben uitgerekend hoe intens het focus 
is (Hoofdstuk 3), wat er gebeurt als de deeltjes in het materiaal bewegen (Hoofdstuk 
7) en hoe groot het focus is (Hoofdstuk 4). Het blijkt dat het focus van, bijvoorbeeld, 
een eierschil precies even scherp is als het focus dat je met de best mogelijke lens zou 
kunnen maken. In Hoofdstuk 6 laten we verder nog zien dat wanordelijke verstrooi- 
ing gebruikt kan worden om meer optische informatie te versturen dan normaal mo- 
gelijk is. 

We zijn ook in staat om licht binnenin een materiaal te focusseren, lets wat zeer 
interessant is voor bijvoorbeeld medische toepassingen. Dit experiment is technisch 
veel lastiger dan licht ergens dwars doorheen focusseren. Om de juiste combinatie 
van fases te vinden is het namelijk nodig om de lichtintensiteit op de gewenste plek 
te meten. Als je, bijvoorbeeld, binnen in de huid wilt focusseren, is het niet mogelijk 
om op die plek een detector aan te brengen. We hebben daarom gebruik gemaakt 
van zeer kleine boUetjes met een fluorescerende verf erin. Dit soort boUetjes zijn wel 
in te brengen. Door nu te meten hoeveel het bolletje fluoresceert, meten we indirect 
hoeveel licht er van buiten af op het bolletje valt. Dit experiment staat beschreven in 
Hoofdstuk 5. 

Verstrooiing van licht lijkt wiskundig gezien erg veel op het transport van elek- 
tronen in een metaaldraad. We kunnen daarom met licht voorspellingen testen die 
gedaan zijn voor elektronische systemen. Het voordeel van onze optische methode 
is dat we het binnenkomende licht heel precies kunnen vormen, lets wat met elek- 
tronen nu niet mogelijk is. Een van de meest verrassende voorspellingen is dat het 
mogelijk moet zijn om een golf te maken die voor de voile 100% wordt doorgelaten 
door een ondoorzichtig materiaal. In Hoofdstuk 8 laten we zien wat dit betekent 
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voor ons experiment: als we de binnenkomende lichtbundel optimaliseren om op 
een punt licht te focusseren moet de intensiteit rondom dat punt ook omhoog gaan. 
Er zal dus in totaal veel meer licht doorgelaten worden (en er zal dus ook veel min- 
der licht worden gereflecteerd). We hebben in een experiment aangetoond dat de 
to tale transmissie inderdaad significant stijgt (Hoofdstuk 9). Hiermee bevestigen we 
een belangrijk theoretisch aspect van bet transport van elektronen, geluid, licht en 
andere golven. 



Conclusie 

We hebben een methode ontwikkeld om met behulp van interferentie laserlicht door 
ondoorzichtige materialen been te sturen. Het op deze manier sturen van laserlicht 
is een nieuw onderwerp in de optica. Dit proefschrift kan daarom het beste gezien 
worden als een eerste verkenning van dit gebied. Het is nu al duidelijk dat er nog 
heel veel verder te experimenteren valt en dat er interessante toepassingen verwacht 
kunnen worden. We hopen dan ook dat dit nieuwe vakgebied ook de komende jaren 
nog tot veel spectaculaire resultaten zal leiden. 
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